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Project  N5ori-07856 
SURVEY 

Several  major  programs  of  calculation  have  been  finished  up  during  the 
period  covered  by  this  Report.  Mr.  Freeman,  in  the  first  place,  has  completed  his 
calculation  on  the  electronic  energy  of  the  OH  molecule,  and  has  fulfilled  all  the 
requirements  for  the  doctor's  degree.  He  leaves  the  group  at  the  end  of  June.  His 
calculations  represent  one  of  the  few  really  complete  treatments  of  a  molecular  wave 
function  by  configuration  interaction;  they  are  particularly  valuable  because  he  has 
carried  through  the  work  both  using  atomic  orbitals,  with  complete  treatment  of  the 
orthogonality  problem,  and  using  molecular  orbitals.  It  is  particularly  interesting 
to  compare  his  results  on  OH  with  the  similar  calculations  of  Scherr  on  N^,  of  Sahni 
on  BH,  and  of  Karo  and  Olson  on  LiH;  these  represent  most  of  the  calculations  of 
comparable  rigor  which  have  been  carried  through  up  to  now.  The  striking  result,  a 
somewhat  disconcerting  one,  is  that  the  calculated  binding  energy  forms  a  very  differ¬ 
ent  fraction  of  the  experimental  value  for  these  various  cases.  One  cannot,  then, 
simply  apply  a  fixed  correction  factor  to  the  calculated  values  to  get  a  correct  value. 
Scherr  and  Freeman  have  both  found  quite  small  fractions  of  the  experimental  value, 
whereas  Sahni  and  Karo  and  Olson  have  found  a  large  fraction  of  the  experimental 
value. 

The  conclusion  which  we  must  draw  from  this,  I  believe,  is  that  further  and 
much  more  extensive  configuration  interaction  must  be  used,  before  one  can  get  cal¬ 
culated  results  that  can  be  trusted.  It  is  well  known  that  though  the  total  energy  of 
the  separated  atoms,  or  of  the  molecule,  can  be  calculated  with  an  error  of  only  about 
half  a  percent  by  the  methods  we  are  now  using,  nevertheless  this  residual  error  in 
absolute  magnitude  is  several  times  as  large  as  the  binding  energy  which  we  are  seek¬ 
ing.  We  should  have  to  get  the  total  energy  with  an  error  of  roughly  a  tenth  as  great 
as  at  present,  before  we  could  begin  to  trust  the  binding  energy  very  seriously.  The 
work  which  Nesbet  and  Watson  have  been  doing  on  configuration  interaction  in  helium, 
the  work  which  I  have  been  doing  on  that  same  problem  (mentioned  in  the  preceding 
Quarterly  Progress  Report),  and  similar  work  which  Shull  and  Lowdin  have  been  doing, 
all  suggest  that  with  a  relatively  simple  configuration  interaction  one  can  decrease  the 
error  of  the  calculation  of  the  helium  atom  to  about  ten  percent  of  the  error  in  the 
Hartree-Fock  energy.  A  similar  improvement  is  required  to  make  the  molecular 
calculations  at  all  reliable,  and  we  must  look  forward  to  similar  configuration  inter¬ 
actions  there,  as  a  goal  to  be  attained  as  rapidly  as  possible. 

The  other  two  calculations  which  have  been  virtually  completed  are  on  energy 
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bands:  Corbato's  work  on  graphite,  and  Howland's  on  KC1.  These  have  both  been 
carried  through  by  the  tight  binding  method,  including  detailed  calculation  of  the  re¬ 
quired  many-center  integrals,  and  both  serve  to  mark  a  new  advance  in  the  use  of  the 
tight -binding  method.  Corbato,  in  his  work  on  the  <r  as  well  as  the  w-bands  in  graph¬ 
ite,  gives  us  a  reliable  estimate  of  the  actual  shape  of  the  occupied  energy  bands  for 
this  important  material.  His  results  could  be  appreciably  in  error  --  he  has  used  an 
assumed  periodic  potential,  without  carrying  through  a  self-consistent  field  calcula¬ 
tion,  and  has  merely  found  the  one-electron  energy  in  this  periodic  potential.  Never¬ 
theless  his  results  are  so  reasonable  that  it  seems  likely  that  more  accurate  calcula¬ 
tions,  when  they  become  possible,  will  merely  verify  the  general  correctness  of  his 
results. 

Howland  presents  in  this  Progress  Report  a  discussion  of  the  lattice  energy 
of  KC1,  but  his  calculations  include  much  additional  material,  as  well,  including  the 
details  of  the  energy  bands.  These  calculations  represent  a  higher  degree  of  sophisti¬ 
cation  than  Corbato's,  in  that  Howland  has  set  up  a  determinantal  wave  function,  and 
computed  its  energy,  including  all  the  various  integrals.  His  values  of  lattice  energy, 
which  can  be  compared  with  Lowdin's  earlier  work,  show  that  Lowdin  omitted  certain 
interactions  which  are  of  significant  size.  These  calculations  provide  one  of  the  most 
accurate  results  which  have  been  obtained  for  the  energy  of  a  crystal  lattice. 

The  reason  why  Howland  was  able  to  get  more  complete  results  than  Corbato 
is  based  essentially  on  the  difference  between  the  two  types  of  crystals.  Graphite  of 
course  is  a  covalently  bonded  crystal,  which  means  that  the  valence  electrons  of 
neighboring  atoms  have  a  very  large  overlap  with  each  other.  Potassium  chloride, 
on  the  contrary,  being  an  ionic  crystal,  has  wave  functions  for  the  outer  electrons 
which  overlap  only  slightly,  their  overlap  taking  care  of  the  repulsion  between  closed 
shells.  It  is  this  much  smaller  magnitude  of  the  overlap  in  KC1  which  makes  it  in¬ 
herently  a  simpler  problem  for  the  application  of  the  tight-binding  method,  and  which 
makes  possible  a  higher  degree  of  approximation  than  in  graphite.  A  calculation  of 
the  order  of  difficulty  which  Corbato  has  carried  through  in  graphite  is  a  very  major 
effort,  which  as  far  as  we  know  has  never  been  even  approached  in  previous  studies 
of  covalently  bonded  crystals  by  the  tight -binding  method.  Corbat6  considers  that 
this  method  would  prove  to  have  very  nearly  insuperable  difficulties  for  a  crystal 
much  more  complicated  than  graphite;  in  this  case  he  was  greatly  helped  by  the  two- 
dimensional  nature  of  the  crystal,  an  almost  unique  property  of  graphite. 

In  addition  to  these  calculations  which  have  reached  virtual  completion, 
Schultz's  considerations  regarding  the  polaron  have  advanced  a  great  deal  in  recent 
months,  and  he  is  bringing  this  work  to  a  conclusion,  hoping  to  finish  it  during  the 
summer.  The  other  work  of  the  group  is  progressing  well,  and  several  other  pieces 
of  work  should  be  brought  to  a  conclusion  before  many  more  months. 


J.  C.  Slater 
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1.  MAGNETIC  HYPERFINE  STRUCTURE  OF  LITHIUM 


A  series  of  calculations  have  been  carried  out  on  the  S  ground  state  of  lithium 

^  This  is  an  attempt  to  obtain  the  mean  value  of  the 


by  methods  described  earlier 
operator 
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The  experimental  value  determined  from  the  magnetic  dipole  hyperfine  structure  is 
<f>  =  2.  90  a"3. 

Configuration  interaction  calculations  have  been  made  using  codetors 
(linear  combinations  of  determinantal  functions)  constructed  from  analytic  orbitals 
orthogonal  to  approximate  self-consistent  orbitals.  Coefficients  in  the  expansion  of 
a  set  of  orthonormal  orbitals  which  includes  the  occupied  self-consistent  orbitals,  con¬ 
structed  from  a  set  of  basic  analytic  orbitals,  have  been  obtained  by  the  Roothaan  pro¬ 
cedure,  ^  under  equivalence  restriction.  This  means  that  the  matrix  equations  of  the 
Roothaan  procedure  are  solved  only  for  orbitals  of  positive  spin;  orbitals  of  negative 
spin  are  constrained  to  have  spatial  factors  identical  with  the  orbitals  of  positive  spin. 
This  introduces  no  approximation  in  a  configuration  interaction  calculation  so  long  as 
all  configurations  which  can  be  constructed  from  the  given  set  of  basic  orbitals  are 
either  included  explicitly  or  are  shown  to  have  a  negligible  effect  on  the  value  of  the 
quantity  being  calculated. 

A  preliminary  calculation  was  made  with  five  s -orbitals,  three  p-orbitals,  and 

two  d-orbitals.  The  calculated  Hartree-Fock  energy  HQO  was  -  7.  43241  a.  u.  and  the 

mean  value  of  the  hyperfine  coupling  operator  f  was  2.  10759  a  3  for  a  single  determinant 

-  3  ^ 

(Hartree-Fock  approximation)  and  2.  66019  a~  after  a  configuration  interaction  calcula¬ 
tion  carried  out  by  second  order  perturbation  methods.  The  contribution  of  configura¬ 
tions  including  p-  and  d-orbitals  to  the  final  value  of  <f>  was  negligible.  It  was  found 
that  configurations  which  were  most  important  in  depressing  the  energy  had  little  effect 
on  <f>  and  vice  versa.  Since  all  coefficients  of  configurations  affecting  <f>  were 
extremely  small  it  does  not  seem  likely  that  significant  errors  can  be  attributed  to 
use  of  the  perturbation  method. 

A  series  of  calculations  were  done  with  six  s -orbitals  and  no  orbitals  of 
higher  angular  momentum,  with  results  as  follows; 

Calculation  Single  determinant  Configuration  interaction 


Hoo 

<f> 

<f> 

2. 

-7.  43250  a.  u. 

2.  05565  a'3 
o 

2. 66706  a 

( 

3. 

-7.43242 

2. 06220 

2.  61702 

4. 

-7.43259 

2. 11881 

2.  67834. 
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{MAGNETIC  HYPERFINE  STRUCTURE  OF  LITHIUM) 

Although  configuration  interaction  gives  a  marked  improvement  over  the 
Hartree-Fock  calculation  of  <f>,  the  results  appear  to  indicate  a  value  of  about 

_3 

2.  68  aQ  for  this  quantity  which  is  insensitive  to  changes  in  the  number  of  basic  or¬ 
bitals,  for  a  wide  range  of  values  of  the  exponent  of  the  sixth  basic  s-orbital.  The 
first  five  orbitals  were  the  same  for  all  calculations.  So  far  it  has  not  been  possible 
to  account  for  the  discrepancy  between  this  result  and  the  experimental  value,  <f>  = 

2. 90  a-3, 
o 

References 

1.  R.  K.  Nesbet,  Quarterly  Progress  Report,  Solid-State  and  Molecular  Theory  Group, 

M.  I.  T.  ,  July  15,  1955,  p.  27. 

2.  R.  K.  Nesbet,  Quarterly  Progress  Report,  Solid-State  and  Molecular  Theory  Group, 

M.  I.  T.  ,  October  15,  1955,  p.  4. 
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2.  ATOMIC  WAVE  FUNCTIONS 


The  most  recent  calculation  for  He  has  given  a  result  of  -2.  9010  a.  u.  This 

involves  an  SCF  calculation  with  4s  like,  3p  like,  2d  like  and  If  like  basis  orbitals,  ^ 

followed  by  diagonalizing  the  20  X  20  configuration  interaction  matrix. 

Whirlwind  has  been  used  in  two  calculations  on  Fe.  Nesbet's^  programs 

6  Z  7  1 

have  again  been  used,  once  on  a  3d  4s  configuration  and  once  on  a  3d  4s  configura¬ 
tion.  The  principal  problem  has  been  that  of  establishing  the  rules  for  choosing  the 
basis  functions.  To  this  purpose,  several  calculations  were  done  on  closed  shell 
transition  elements  for  which  there  are  Hartree-Fock  results.  Good  agreement  with 
the  H-F  results  was  the  principal  criterion  for  the  basis  orbitals.  A  good  choice  of 
basis  orbitals  is  very  important,  since  the  size  of  the  computer  limits  the  number 
which  may  be  used.  For  example,  we  are  using  three  basis  functions  for  the  Fe  3d 
orbitals.  It  appears  that  the  best  choice  of  basis  functions  in  the  closed  shell  case 
does  not  exactly  interpolate  to  a  best  choice  for  the  unfilled  shell  case  and  we  are 
still  in  the  process  of  forming  the  rules  for  these  choices. 

Reference 

1.  R.  K.  Nesbet,  Quarterly  Progress  Report,  Solid-State  and  Molecular  Theory  Group, 
M.I.  T.  ,  April  15,  1955,  p.  38. 
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3.  THE  LITHIUM  HYDRIDE  MOLECULE 


Calculations  of  the  electronic  energy  for  the  ground  state  of  lithium  hydride 
by  the  self-consistent  molecular  orbital  method  are  in  progress.  The  electronic  con¬ 
figuration  of  this  state  is  represented  as  (lor)2  (2<r)2  T.  •  The  binding  energy  has  been 
obtained  for  five  internuclear  distances,  and  the  results  are  given  in  Table  3-1. 

Table  3-1 

Internuclear  distance  (atomic  units) 


2.  0 

2.  6 

3.  0 

3.  5 

4.  0 

SCF-MO  Method: 

+  1.  1413 

-0. 8702 

-1.  1989 

-1.  1464 

-0. 8371 

Case  1: 

+  1.  1019 

-0. 9310 

-1.  2794 

-1.  2912 

-1.  1072 

A: 

+0.  0394 

+0. 0608 

+0.  0805 

+0. 1448 

+  0.  2701 

The  molecular  orbitals  are  found  as  a  linear  combination  of  the  four  atomic 
orbitals  which  we  are  using.  An  inspection  of  the  coefficients  of  the  atomic  orbitals 
forming  the  molecular  orbitals  for  internuclear  distances  near  the  equilibrium  value 
indicates  that  the  l<r  corresponds  to  a  nearly  pure  lithium  Is  orbital,  and  that  the  Z<r 
corresponds  on  the  whole  to  a  combination  of  the  lithium  2s  and  hydrogen  Is  orbitals, 
with  a  lesser  mixing  of  the  lithium  2p  function.  Therefore,  it  appears  that  an  approxi¬ 
mate  comparison  can  be  made  with  prior  work  (case  1  of  the  previous  Progress  Report) 
in  which  we  excluded  2p  orbitals  in  the  configuration  interaction  problem.  In  Table 

we  have  given  these  results  again  and  also  the  difference  in  the  binding  energy  cal¬ 
culated  by  the  two  methods. 

We  may  note  that  the  single  configuration  molecular  orbital  treatment  is  be¬ 
coming  worse  as  the  internuclear  distance  increases.  Although  we  have  only  made  a 
very  preliminary  analysis  of  the  effects  of  configuration  interaction,  we  should  expect 
that  the  final  result  for  the  binding  energy  will  be  about  that  obtained  for  case  2  of 
the  preceding  Progress  Report,  and  perhaps  somewhat  better  since  we  have  the  addi¬ 
tional  mixing  in  of  the  Is  lithium  function.  Our  previous  work,  however,  has  indicated 
that  the  inclusion  of  the  Is  lithium  electrons  will  be  comparatively  unimportant. 

The  four  molecular  orbitals  which  are  constructed  in  the  process  of  making 
the  calculation  are  of  <r  symmetry  (zero  component  of  orbital  angular  momentum  about 
the  internuclear  axis),  since  we  have  not  as  yet  introduced  the  2p(m  =  +  1)  or  the 
2p(m  =  -  1)  atomic  orbitals  into  the  calculation.  It  now  seems  that  it  would  be  of  suffi¬ 
cient  interest  to  include  these  2p  functions  in  both  of  the  schemes  which  we  are  using. 

In  the  Lowdin  formalism  this  would  include  an  additional  configuration  to  the  six  pre¬ 
viously  used  and  discussed  as  case  2.  ^  The  Roothaan  method  goes  through  straight - 
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forwardly  with  the  proper  listing  of  the  additional  integrals  in  the  Whirlwind  programs. 
The  radial  part  of  the  integrations  have  already  been  performed  so  that  only  the  angu¬ 
lar  factors  need  to  be  computed. 

Computation  of  the  integrals  for  the  lithium  hydride  molecule  at  an  internu- 
clear  distance  of  eight  atomic  units  is  very  nearly  completed.  We  have  found  that  the 
convergence  of  hybrid  and  exchange  integrals  is  slow,  but  still  satisfactory  even  at 
this  distance.  The  binding  energy  will  be  computed  for  this  distance  using  the  con¬ 
ventional  configuration  interaction  scheme  as  well  as  the  self-consistent  molecular 
orbital  procedure.  Comparison  of  these  results  should  prove  interesting  because  of 
the  failure  of  the  ground  state  configuration  molecular  orbital  treatment  at  larger 
internuclear  distances. 
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4.  ELECTRONIC  ENERGY  OF  THE  OH  MOLECULE 


This  is  the  final  report  under  this  title.  The  calculations  described  only 
briefly  in  previous  reports  have  been  completed  and  submitted  as  a  thesis  to  the 
Physics  Department.  In  conclusion,  some  of  the  more  pertinent  details  omitted  in 
previous  contributions  to  these  reports  along  with  some  of  the  results  of  the  configura¬ 
tion  interaction  (Cl)  studies  are  presented  here. 

Cl  Calculation  Using  Non-Orthogonal  AO's 

Under  this  heading,  we  have  described  in  previous  reports,  a  Cl  study  which 
used  as  basis  those  orbitals  natural  to  the  isolated  atoms  or  AO's.  The  non-orthogonality 
of  these  AO's  and  of  the  various  configurations  formed  from  them  was  treated  without 
approximation  by  the  method  of  Lowdin.  The  necessary  integrals  involving  the  AO's 
were  calculated  using  methods  described  previously.  Table  4-1  lists  the  AO  configura¬ 
tions  of  2ir  symmetry  enumerated  in  a  previous  Report.  ^  In  this  table  s  is  the  oxygen 
Is  function,  S  the  oxygen  2s  function,  and  pQ,  p+1  p  denote  the  three  2p  orbitals  in 
oxygen  with  components  of  angular  momentum  equal  to  respectively  0,  +1,  and  -1. 
h  denotes  the  hydrogen  Is  state  which  is  not  made  orthogonal  to  the  oxygen  orbitals. 

Table  4-1 

Enumeration  of  the  AO  configurations 


State 

Coefficient 

s 

s 

h 

Po 

p+ 

P. 

*1 

1 

+  - 

+- 

+  - 

+  - 

+ 

*2 
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+- 

+- 
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+  - 
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+ 
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i/SI 
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(+)  and  (-)  denote  the  ordinary  spin  functions.  Each  row  represents  a  single  Slater 
determinant,  a  state  being  taken  as  a  sum  of  rows  with  coefficients  (given  to  the  left 
of  each  row)  which  make  these  states  eigenfunctions  of  the  total  spin  operator. 

In  Table  4-2,  the  expansion  coefficients,  C^,  and  the  lowest  root  of  the 
secular  equation  giving  the  binding  energy,  Eg,  are  listed  both  for  the  five  configura¬ 
tions  formed  by  keeping  the  Is  and  2s  oxygen  orbitals  doubly  filled,  and  for  the  more 
complete  ten  configuration  treatment  which  results  when  the  2s  electrons  are  excited. 
From  Table  4-2  we  observe  that  the  three  most  important  states  are  and  Jg; 

only  corresponds  to  a  polar  state  OH+. 

None  of  the  configurations  have  a  diagonal  energy  lower  than  the  total  energy, 

E  ,  of  the  separated  atoms;  in  particular,  this  is  also  true  for  the  configurations 
and  which  correspond  to  valence  bond  states  (since  they  describe  bonds  between 
the  pQ  and  h  and  the  S  and  h  orbitals,  respectively).  A  linear  combination  of  these 
states,  forming  an  s-pQ  hybridized  bond  with  the  h  function,  still  fails  to  give  binding 
even  though  the  coefficients  are  determined  by  a  variational  procedure.  This  is  to  be 
expected  since  the  atoms  in  this  hydride  are  fairly  close  together,  and  so  are  in  that 
region  in  which  the  familiar  Heitler-London  approximation  is  least  valid.  Hence  a 
description  in  terms  of  rigid  AO's  such  as  the  ones  employed  here  should  be  a  poor 
approximation. 

Eg  =  -.  0627  Rydbergs  for  the  5x5  matrix  of  interaction;  the  remaining 
five  states  succeed  in  lowering  it  to  Eg  =  -.  0764  Rydbergs.  Since  the  observed  value 
is  -.  337  Rydbergs  (including  the  zero-point  energy  of  .017  Rydbergs)  the  ten  configura¬ 
tions  gives  only  22  percent  of  the  experimental  value.  (More  will  be  said  later  about 
this  result,  when  the  MO  scheme  has  been  presented.  )  The  total  energy  is  of  course 
in  better  percentage  agreement  with  experiment  yielding  99.  5  percent  of  the  measured 
value.  The  lowering  of  the  energy,  brought  about  by  those  states  formed  by  exciting 
the  2s  electrons,  is  very  small  --  hardly  justifying  the  amount  of  computational 
effort  that  was  involved. 

Table  4-2 

Values  of  the  expansion  coefficients 


States  ijij  through 


States  through 


-  . 0627  (Ryd.  ) 

1 

II 

to 

w 

.  0764  (Ryd. ) 

=  . 011209 

Cj  =  . 014767 

c6  =  - 

005303 

=  . 016577 

c ^  =  .  017306 

c  7  =  "■ 

010048 

=  . 471287 

c3  = .414435 

II 

00 

a 

122771 

=  . 626745 

c4  =  . 607075 

c9  =  - 

02605 

=  . 036056 

c5  =  . 035611 

o 

*— • 

°« 

1 

02723 
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Since  this  calculation  is  the  first  of  its  kind,  comparison  can  be  made  only 
with  the  results  obtained  by  other  theoretical  methods.  In  particular,  since  we  have 
worked  out  the  same  problem  by  the  MO-LCAO  method,  we  have  in  essence  built  in 
our  own  check  on  the  results  just  described. 

Cl  Study  Using  MO  LCAO's 

In  this  scheme  another  configuration  interaction  calculation  is  carried  out 

using  as  basis  a  set  of  orthonormal  MO's  which  are  determined  from  the  same  fixed 

TZT — 

set  of  AO's  used  in  the  last  section.  These  MO's  are  determined  by  Roothaan's  ' 

SCF  procedure  for  finding  the  "best"  set  of  MO  LCAO's  which  minimize  the  energy 
of  a  single  configuration  at  each  value  of  the  internuclear  separation.  In  this  sense 
the  MO's  are  not  of  the  rigid  form  employed  in  the  AO  scheme.  These  SCF  MO's 
will  then  be  used  as  basis  functions  for  a  configuration  interaction  calculation  of  the 
molecular  ground  state.  In  contradistinction  to  the  work  presented  in  the  previous 
section  no  elaborate  formalism  for  handling  the  matrix  elements  is  needed.  The  MO's 
form  an  orthonormal  set  and  so  none  of  the  orthogonality  difficulties  encountered 
there  appear  in  this  study. 

Determination  of  the  "Best"  Single  Configuration  MO's 

The  LCAO  approximation  constructs  symmetry  orbitals  from  a  set  of  atomic 
functions.  If  we  call  the  AO's  <j>j,  then  the  MO's,  are  chosen  in  the  form 

=  I  Cy  4>.  (4-1) 

The  C--  are  expansion  coefficients,  undetermined  as  yet.  To  satisfy  the  Pauli  prin- 
ciple,  we  shall  approximate  the  molecular  state  function  by  a  single  Slater  determinant, 
made  up  from  the  as  MO's.  We  want  to  determine  that  set  of  MO's  which  will 
give  the  lowest  possible  energy  for  subject  to  the  restriction  that  the  tpAs  form  an 
orthonormal  set.  These  "best"  MO's  are  determined  by  the  Roothaan  method  which 
uses  the  variational  procedure  in  a  self-consistent  way  to  select  the  coefficients,  Cy, 
which  make  stationary.  Since  Meckler  and  Nesbet  have  discussed  their  programs 
for  Whirlwind  for  carrying  out  these  computations  in  other  issues  of  these  Reports, 
we  shall  give  no  further  details  here. 

For  OH,  from  the  six  available  atomic  functions  only  the  following  twelve 
SCF  LCAO  MO's  can  be  formed:  l<r,  2a,  3a,  4a,  w  ,  w’  of  plus  spin  and  Ta,  2a,  3a, 

?a,  ir+,  w’,  of  minus  spin.  These  symmetry  MO's,  as  indicated  by  the  notation,  be¬ 
long  to  irreducible  representations  of  the  group  C  ,  the  symmetry  group  of  the 

oov 

molecule;  the  orbitals  are  labeled  in  the  order  of  increasing  energy.  Since  the  OH 
molecule  has  an  odd  number  of  electrons  in  a  non-closed  shell  ground  state,  the  MO's 
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of  different  spin  cannot  be  identical,  and  so  must  have  different  LCAO  coefficients. 

This  is  due  to  the  difference  in  exchange  energy  between  those  electrons  having  the 
same  spin  as  the  odd  electron  and  those  of  opposite  spin.  The  true  solution  of  the 
H-F  equations  would  therefore  result  in  different  sets  of  integro-differential  equa¬ 
tions  for  symmetry  orbitals  of  opposite  spin.  Furthermore,  the  standard  eigenvalue 
form  of  the  H-F  equations  could  no  longer  be  retained.  For  these  reasons  we  chose 
not  to  solve  the  more  exact  (and  difficult)  H-F  equations  for  the  MO's  of  both  spins 
but  arbitrarily  define  the  LCAO  MO's  to  be  identical  for  electrons  of  either  spin  and 
solve  the  SCF  equations  for  one  spin  set  only.  These  approximate  solutions  of  the 
(now  modified)  Roothaan  procedure  satisfy  the  equivalence  restriction  of  Nesbet  and 
are  thus  insured  of  forming  an  orthonormal  set  of  symmetry  orbitals  with  which  to 
approximate  the  molecular  ground  state.  While  not  exactly  equivalent  to  the  true  solu¬ 
tions  of  the  H-F  equations  (which  we  cannot  find),  our  MO's  should  differ  from  these 
only  in  a  negligible  way. 

The  ground  state  configuration  of  OH  (with  2ir  symmetry)  is  taken  to  be 

0<r)2  (2«r)2  (3a)2  (/)2  (*' a),  (4-2) 

in  which  each  of  the  MO's  is  doubly  occupied  except  for  the  ir  orbital  which  alone  has 
an  unpaired  (a)  spin.  The  4a  orbitals  are  not  occupied  in  the  ground  state  and  so 
represent  excited  orbitals.  The  MO's  expressed  as  LCAO's  are  simply  written  as 


ma 

+ 

ir 


C  s  +  C  S  +  C  p+  C„,  h;  m'1,2,3,4 

irij  m3ro  m4 


TT 


=  P 


(4-3) 


The  degenerate  orbitals,  ir+  and  ir',  consisting  of  only  one-member  terms,  have 
their  LCAO  coefficient,  uniquely  determined  by  the  normalization  condition.  The  H-F 
equations  thus  reduce  to  a  standard  eigenvalue  problem  involving  the  solution  of  a 
4X4  secular  determinant  for  the  LCAO  MO's.  Table  4-3  lists  the  computed  coeffi¬ 
cients  for  the  SCF  LCAO's  along  with  the  values  of  the  orbital  energies,  e  ,  for 
three  values  of  the  internuclear  distance.  The  inner-outer-shell  mixing  of  the  AO's 
is  very  small.  The  la  MO  is  made  up  almost  purely  from  the  s  orbital  (the  oxygen 
atomic  Is  function),  the  higher  a  orbitals  contain  only  a  small  contribution  from  this 
AO.  The  variation  in  charge  distribution  among  the  AO’s  is  shown  in  Table  4-3  both 
as  a  function  of  internuclear  distance  and  among  the  a  orbitals  for  a  particular  inter¬ 
nuclear  distance. 
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Note:  These  one-electron  energies  are  for  the  MO's  with  plus  spin  only.  The  MO's 
with  negative  spin  have  energies  which  differ  from  these  by  an  exchange  inte¬ 
gral  with  the  it"  electron.  The  unit  of  energy  is  2  Rydbergs. 
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The  total  energy  of  the  molecule  is  the  sum  of  the  nuclear  repulsion  energy 
and  the  one-electron  energies  of  the  occupied  SCF  LCAO's  in  the  configuration. 

The  binding  energy  is  calculated  by  subtracting  the  total  energy  of  the  free  atoms 
from  the  value  computed  for  the  molecule.  In  Table  4-4,  the  results  of  these  calcu¬ 
lations  are  presented  and  compared  with  the  corresponding  experimental  values.  The 
total  molecular  energy  is  quite  good  whereas  the  binding  energy  seems  a  bit  disap¬ 
pointing.  Actually  though,  these  values  are  consistent  with  the  results  of  other  SCF 

(3) 

LCAO  calculations.  For  comparison,  we  list  in  Table  4-5  the  results  of  Sahni  on 
(4) 

BH  and  Seherr'  '  on  N^,  the  only  other  molecules  reported  to  have  been  computed 
employing  no  approximations  other  than  those  inherent  in  the  method.  (The  results 
of  Karo  and  Olson  on  LiH  discussed  by  them  in  this  Report  were  not  available  when 
Table  4-5  was  prepared  and  so  could  not  be  included  therein.  )  Both  authors  used 
Slater  AO's  with  Slater's  values  for  the  exponents;  the  energy  was  calculated  at  a 
value  of  the  internuclear  distance  equal  to  the  experimentally  observed  equilibrium 
value.  The  data  for  OH  listed  in  Table  4-5  was  calculated  at  the  internuclear  separa¬ 
tion  equal  to  2.  0.  The  total  energy  of  the  separated  atoms  for  OH  is  seen  to  be  better 
than  for  N^  or  BH.  This  is  to  be  expected,  since  Hartree-Fock  functions  have  more 
physical  content  as  AO's  than  do  the  simple  Slater  functions.  Correspondingly,  the 
molecular  energy  also  compares  just  as  favorably.  The  magnitude  of  the  binding 

energy  is  about  the  same  for  all  three  calculations.  Sahni  uses  the  value  of  2.  6  e.  v. 

(5) 

for  the  observed  dissociation  energy.  Since  Gaydon  '  actually  recommends  the  value 
of  3.  0  t  .  4  e.  v.  ,  both  values  were  included  in  Table  4-5  to  indicate  the  range  of  the 
probable  binding  energy  which  was  obtained  for  BH.  Sahni  does  not  specifically  say 
that  he  included  a  correction  for  the  zero  point  vibration  energy  and  so  it  is  probably 
quite  correct  to  assume  that  no  such  correction  was  made.  If  included,  this  would 
decrease  his  binding  energy  result,  bringing  it  more  in  line  with  the  other  results 
quoted  in  Table  4- 5b. 

As  a  check  on  the  stability  of  the  binding  energy  with  respect  to  changes  in 
values  of  the  integrals,  the  oxygen  kinetic  energy  integrals  were  changed  by  as  much 
as  .  04  Rydbergs.  (These  integrals  were  chosen  because  they  represent  the  most  in¬ 
accurate  integrals  calculated  in  this  work. )  The  total  energies  of  both  the  molecular 
and  free  atom  systems  changed,  but  the  binding  energy  was  almost  completely  un¬ 
affected. 

Cl  Using  SCF  LCAO's 

Having  determined  that  set  of  LCAO  MO's  which  gave  the  lowest  energy  of 
a  single  configuration,  we  are  now  in  a  position  to  study  the  effects  of  configuration 
interaction  in  further  lowering  the  energy  of  the  molecule.  All  possible  configurations 
possessing  the  appropriate  symmetry  (^w)  are  to  be  formed  from  the  available 
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LCAO  MO's  in  Eq.  (4-3).  Since  only  the  4<r  orbitals  were  not  included  in  the  con¬ 
figuration,  the  excited  configurations,  are  easily  formed  by  promoting  electrons 
from  the  lower  lying  occupied  orbitals  into  the  excited  4<r  states.  These  are  presented 
in  Table  4-6  with  the  same  notation  as  was  used  in  Table  4-1.  The  ^  are  orthonormal 
since  the  orbitals  are,  and  so  the  formulae  for  the  matrix  elements  of  the  Hamiltonian 
between  two  Slater  determinants  are  much  simpler  than  for  a  non-orthogonal  basis. 

The  Cl  problem  is  solved  by  diagonalizing  the  energy  matrix  (H^^)-  Since 
there  are  17  states  made  up  from  23  Slater  determinants,  the  computation  of  all  the 
matrix  elements  between  them  is  no  small  job.  Fortunately,  the  simplifications  cited 
by  Nesbet^  were  used  to  reduce  substantially  the  magnitude  of  the  calculations. 

According  to  Brillouin's  theorem,  the  matrix  elements  of  energy  between  de¬ 
terminants  which  differ  by  a  single  orbital  are  identically  zero.  Our  set  of  LC/tO 
MO's  are  solutions  of  the  H-F  equations,  modified  by  the  equivalence  restriction. 

That  is,  for  single  substitutions  of  orbitals  of  plus  spin  Brillouin's  theorem  is  auto¬ 
matically  satisfied  whereas  the  single  substitution  matrix  elements  involving  orbitals 
of  negative  spin  are  not  zero,  but  are  found  nevertheless  to  be  very  small.  We  need 
therefore  only  include  those  configurations  which  are  formed  by  double  substitution 
for  those  occupied  orbitals  in  This  immediately  eliminates  the  states  ij>7  through 
$  j0  inclusive  from  the  configuration  matrix. 

The  contributions  of  the  excited  configurations,  towards  lowering  the 
energy  of  $0>  can  at  best  be  very  small,  since  as  we  have  seen,  itself  is  a  fairly 
good  representation  of  the  molecular  ground  state  near  the  equilibrium  internuclear 
distance.  (Cl  effects  become  increasingly  important  with  increasing  internuclear 
distance  because  the  MO's  become  an  every  poorer  approximation  as  the  molecular 
system  approaches  the  free  atom  stage. )  For  this  reason  Nesbet^  has  proposed  the 
use  of  second  order  perturbation  theory  for  approximately  diagonalizing  the  energy 
matrix.  Using  the  familiar  formulae  of  perturbation  theory,  the  configurations 
and  tjig  through  J ^  must  play  a  negligible  role,  since  the  1  <r  state  lies  about  40  Ryd¬ 
bergs  below  the  other  one-electron  states. 

If  the  second  order  perturbation  formula  is  used  to  diagonalize  the  secular 

equation  the  total  energy  (and  hence  the  dissociation  energy  as  well)  is  lowered  by 

.  0323  Rydbergs  at  the  internuclear  distance  equal  to  1.  8  (which  is  the  value  used  in 

the  Cl  AO  method  described  previously  in  the  first  section);  the  bordered  determinant 
♦ 

approximation  results  in  a  similar  lowering  of  .  0  320;  the  exact  diagonalization  of 
the  6x6  secular  equation  on  Whirlwind  gives  a  value  of  .  0303  for  the  same  quantity. 
Both  the  perturbation  theory  method  and  the  bordered  determinant  method  overestimate 


xhis  method  diagonalizes  the  matrix  as  approximated  by  the  diagonal  plus  first 
column  matrix  elements. 


(■ 


& 
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Table  4-6 

Enumeration  of  the  MO  configurations 


State 

Coefficient 

lor 

2o- 

3<r 

+ 

IT 

tr 

4a- 

*o 

1 

+  - 

+  - 

+  - 

+  - 

+ 

*1 

1 

+  - 

+  - 

+  - 

+ 

+  - 

*2 

1 

+  - 

+  - 

+  - 

+ 

+  - 

*3 

1 

+  - 

+  - 

+  - 

+ 

+  - 

t  a 
*4 

1 

+  - 

+  - 

+  - 

+ 

+  - 

h 

i  M* 

+  - 

+ 

- 

+  - 

+ 

+  - 

-1/-/2 

+  - 

- 

+ 

+  - 

+ 

+  - 

h 

2/76 

+  - 

+ 

+ 

+  - 

- 

+  - 

-  1/76 

+  - 

+ 

- 

+  - 

+ 

+  - 

-  1/76 

+  - 

- 

+ 

+  - 

+ 

+  - 

T  b 

h 

1/72 

+  - 

+  - 

+ 

+  - 

+ 

- 

-i/TT 

+  - 

+  - 

- 

+  - 

+ 

+ 

T  b 

2/76 

+  - 

+  - 

+ 

+  - 

- 

+ 

-  1/76 

+  - 

+  - 

+ 

+  - 

+ 

- 

-  1/76 

+  - 

+  - 

- 

+  - 

+ 

+ 

T  b 
*9 

1/72 

+  - 

+ 

+  - 

+  - 

+ 

- 

-  772 

+  - 

- 

+  - 

+  - 

+ 

+ 
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Table  4-6  (con'd) 

Enumeration  of  the  MO  configurations 


State 

Coefficient 

1  (T 

2(r 

3<r 

+ 

TT 

TT 

4o- 

T  b 
*10 

2/  y? 

+  - 

+ 

+  - 

+  - 

- 

+ 

-  l/7? 

+  - 

+ 

+  - 

+  - 

+ 

- 

-  l/7? 

+  - 

- 

+  - 

+  - 

+ 

+ 

T  a 

in 

i/yi 

+ 

+  - 

+  - 

+  - 

+ 

- 

-  i/7T 

- 

+  - 

+  - 

+  - 

+ 

+ 

i  a 
*12 

2/  y? 

+ 

+  - 

+  - 

+  - 

- 

+ 

-  i/yr 

+ 

+  - 

+  - 

+  - 

+ 

- 

-  1/7? 

- 

+  - 

+  - 

+  - 

+ 

+ 

ilB3 

i/yr 

+ 

+  - 

- 

+  - 

+ 

+  - 

-  l/TI 

- 

+  - 

+ 

+  - 

+ 

+  - 

*14" 

2/76 

+ 

+  - 

+ 

+  - 

- 

+  - 

- 1  /  y? 

+ 

+  - 

- 

+  - 

+ 

+  - 

- 1/7? 

- 

+  - 

+ 

+  - 

+ 

+  - 

ii5a 

i/7T 

+ 

- 

+  - 

+  - 

+ 

+  - 

-  1/7? 

- 

+ 

+  - 

+  - 

+ 

+  - 

ii6a 

2/7? 

+ 

+ 

+  - 

+  - 

- 

+  - 

- 1/7? 

+ 

- 

+  - 

+  - 

+ 

+  - 

- 1/7? 

- 

+ 

+  - 

+  - 

+ 

+  - 

.  Negligible  by  second-order  perturbation  theory. 
.  Negligible  by  Brillouin's  theorem. 
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the  Cl  and  so  where  possible  exact  diagonalization  of  the  secular  equation  is  pre¬ 
ferred. 

The  result  then  of  solving  the  Cl  problem  in  terms  of  the  MO  LCAO's  as 
basis  is  to  give  a  value  of  the  binding  energy  equal  to  092  Rydbergs.  This  compares 
well  with  the  value  of  -.  0764  Rydbergs  reported  for  the  AO  scheme  in  the  first  sec¬ 
tion.  This  comparison  can  only  be  made  at  the  nuclear  separation,  R,  equal  to  1.  8 
at  which  the  AO  calculation  was  carried  out.  The  minimum  of  energy  for  the  single 
determinant  approximation  (see  Table  4-4)  occurs  at  R  =  2.  05  atomic  units  which  is 
some  ten  percent  larger  than  the  observed  value  of  1.  834  a.  u.  This  result  indicates 
that  had  the  Cl  calculation  of  the  first  section  been  carried  out  at  a  value  of  nuclear 
separation  larger  than  the  value  actually  used,  we  would  certainly  have  obtained  a 
lower  binding  energy  than  was  actually  found. 

Since  the  set  of  MO's  is  formed  from  the  basic  set  of  atomic  orbitals,  the 

same  subspace  of  Hilbert  space  is  spanned  by  both  sets.  If  all  possible  configurations 

were  formed  from  either  set,  the  n-electron  system  would  be  equivalently  described 
(7) 

by  both.  '  Neither  would  form  a  more  complete  description  than  the  other  and  so  the 
same  energy  levels  would  be  derived  from  each. 

The  two  calculations  described  cannot  give  the  same  results  because  the  Cl 
expansions  are  not  equivalent.  The  configurations  formed  from  the  AO's  as  basis 
(see  Table  4-1)  did  not  include  those  states  which  resulted  from  exciting  the  oxygen 
Is  electrons  into  higher  lying  states.  The  Cl  calculation  involving  the  SCF  LCAO's, 
expanded  the  total  wave  function  as  a  linear  combination  of  all  the  configurations  which 
could  be  obtained  from  the  basic  set  of  MO's.  True,  the  configurations  resulting 
from  exciting  la  electrons  into  higher  orbitals  were  found  to  come  into  the  expansion 
with  extremely  small  coefficients  and  so  were  completely  left  out  of  the  secular  de¬ 
terminant.  But  the  other  a  orbitals  each  contain  a  small  amount  of  the  Is  oxygen  AO, 
and  so  the  configurations  already  contain  a  small  but  appreciable  contribution  from 
the  Is  states.  If  there  were  no  mixing  in  the  LCAO's  between  the  s  orbital  and  the 
other  spherically  symmetric  AO's,  or  if  all  the  AO  configurations  were  used,  the  two 
solutions  would  then  agree.  The  value  of  -.  0764  Rydbergs  obtained  for  the  binding 
energy  by  the  AO  Cl  scheme  is  only  .  016  Rydbergs  above  the  result  of  the  SCF  LCAO 
configuration  interaction  calculation.  This  seems  to  be  a  reasonable  amount  to  attrib¬ 
ute  to  the  Cl  lowering  which  would  be  introduced  by  those  (seven)  configurations 
neglected  in  the  AO  scheme.  It  may  therefore  be  concluded  that  the  results  are  con¬ 
sistent  among  themselves  and  seem  to  indicate  a  fair  degree  of  numerical  accuracy. 

Conclusions 

The  results  of  the  Cl  studies  are  not  as  good  as  we  would  like  them  to  be. 

For  closer  agreement  with  experiment  we  should  set  up  a  more  unrestricted  Cl 
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calculation  by  expanding  the  basic  set  of  functions,  mixing  in  (say)  AO's  with  higher 

(8) 

orbital  angular  momentum  to  provide  the  possibility  of  increased  angular  correlation.  '  ' 
Such  an  expanded  Cl  calculation  would  completely  alter  the  character  of  the  problem  we 
have  completed  here.  There  are  no  other  Hartree-Fock  orbitals  for  atomic  oxygen 
available  with  which  to  augment  the  set  used  here.  If  more  orbitals  are  to  be  introduced, 
it  would  be  best  to  choose  an  entirely  new  set  for  the  atom,  so  chosen  as  to  provide 
the  maximum  degree  of  correlation. 

The  major  inadequacy  of  the  calculation  lies  in  what  seems  to  be  a  poor  de¬ 
scription  of  the  n  electrons.  This  is  most  easily  seen  in  the  MO  picture  where  the  ir 
electrons  are  described  by  one-electron  functions  which  are  identical  to  the  p+  and  p 
oxygen  orbitals.  In  the  molecular  environment  this  picture  of  the  outer  electrons  can¬ 
not  be  right.  The  complete  localization  of  the  n  electrons  around  the  oxygen  center  is 
correct  only  at  infinite  separation.  In  the  molecular  environment  some  of  the  ir-electron 
charge  clouds  must  be  described  by  functions  which  spend  an  appreciable  time  near 
the  hydrogen  center.  Our  description  cannot  do  this,  since  we  are  using  a  restricted 
set  of  AO's.  We  are  again  forced  to  conclude  that  in  order  to  improve  our  results  we 
must  expand  our  basic  set  as  outlined  above. 
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5.  GROUND  STATE  OF  THE  HYDROGEN  MOLECULE 


The  exchange  and  hybrid  integrals  which  were  evaluated  by  Merryman's  pro¬ 
gram  were  found  to  be  inaccurate  for  interatomic  separation  of  one  and  two  atomic 
units.  These  integrals  are  now  being  recalculated  by  the  method  of  expansion  about 
another  center  which  is  programmed  for  Whirlwind  by  F.  J.  Corbato.  ^  Expansion 
of  the  exchange  integrals  gives  infinite  series  which  converge  rapidly  for  all  the 
required  integrals.  For  the  internuclear  distance  of  one  atomic  unit  it  was  found  that, 
in  all  cases,  seven  terms  of  the  series  gave  accuracy  of  six  figures.  The  combination 
of  terms  in  the  series  was  done  on  Whirlwind  by  a  program  written  by  L.  P.  Howland. 
The  Whirlwind  time  for  evaluation  of  the  final  integrals  averaged  about  70  seconds  per 
integral. 
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6.  POLARIZATION  EFFECTS  IN  THE  FLUORINE  ION 


The  numerical  integration  of  the  equations  for  our  perturbed  functions 
still  show  considerable  sensitivity  to  the  smoothness  of  the  potential  and  the  mesh 
size.  However,  for  the  uniform  field  case  we  have  obtained  a  complete  check  on  the 
perturbed  functions  of  d-like  symmetry  and  considerable  improvement  in  the  p  and 
s-like  ones.  Also,  we  have  recently  obtained  the  v  function  of  d  symmetry  for  the 
dipole  term  of  a  point  charge  perturbing  potential  at  a  finite  distance,  R,  from  the 
origin.  For  R  =  1.  76  atomic  units  the  v  function  for  Neon  has  a  behavior  near  the 
origin  very  similar  to  the  corresponding  uniform  field  (R  =  ao)  case  and  a  maximum  at 
approximately  r  =  R. 

Because  of  the  delays  and  difficulties  in  obtaining  our  distorted  functions  by 
direct  integration  we  have  undertaken  a  calculation  for  the  HF  molecule  based  on  a 
small  set  of  unperturbed  functions  to  describe  the  F  ion  and  including  a  single  d-like 
distorted  function.  This  d  function  has  been  chosen  to  have  a  maximum  very  close  to 
the  position  of  the  H+  point  charge.  For  the  unperturbed  functions  the  Is  is  F.  W. 
Brown's  fit  of  his  Hartree-Fock  F  calculation,  the  2s  is  the  single  Slater  atomic  orbital 
which  most  nearly  fits  Brown's  2s  function  when  orthogonalized  to  the  Is,  and  a  single 
Slater  AO  was  used  for  the  2p  with  the  exponent  of  the  exponential  chosen  so  the  AO 

has  the  same  maximum  as  the  Hartree-Fock  function.  Explicitly  the  radial  parts  of 

.  .  ,  ,  ..  -8.  75r  -2. 40r  -2. 67r  ,  2  -2.  OOr 

the  Is,  2s,  2p,  d-distorted  functions  are:  e  ,  re  ,  re  ,  and  r  e  , 

respectively. 

We  have  obtained  all  of  the  necessary  atomic  and  potential  integrals  from  R.  K. 
Nesbet's  Whirlwind  program  (described  by  the  writer  and  R.  K.  Nesbet  in  previous  Prog¬ 
ress  Reports).  These  have  been  combined  with  the  appropriate  coefficients  and  classified 
in  terms  of  the  molecular  symmetry.  That  is,  the  axis  of  quantization  for  the  atomic  func¬ 
tions  is  chosen  as  the  internuclear  axis  allowing  the  atomic  quantum  number  m^to  go  directly 
into  the  corresponding  molecular  quantum  number  with  the  combining  coefficients 
the  C  (im^,  i'm^'J's  tabulated  by  Condon  and  Shortley.  v  '  Using  these  integrals  the 
parameter  tapes  for  the  Whirlwind  program  of  the  Roothaan  scheme  are  being  prepared. 

Reference 

1.  E.  U.  Condon  and  G.  H.  Shortley,  The  Theory  of  Atomic  Spectra,  Cambridge 

University  Press,  1953,  p.  178. 
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7.  ENERGY  BANDS  OF  GRAPHITE 


The  computation  of  energy  bands  in  a  two-dimensional  model  of  graphite  by 
means  of  the  tight-binding  method  has  been  completed.  Current  effort  is  being  devoted 
to  the  thesis  write-up  of  this  work.  Thus,  the  following  report  is  only  a  brief  resume 
of  the  basis  of  the  calculation  and  the  major  results  which  were  obtained. 

The  graphite  energy  band  calculations  were  carried  out  assuming  a  structure 
of  an  infinite  two-dimensional  hexagonal  lattice  of  Carbon  atoms.  The  atomic  orbitals 
used  to  make  up  the  one-electron  Bloch  waves  were  the  Hartree-Fock  atomic  orbitals 
of  the  ground  state  neutral  Carbon  3P  configuration  calculated  by  Jucys.  ^  The  Is 
and  2p  orbitals  were  each  fitted  by  a  linear  combination  of  three  corresponding  analytic 
Slater  atomic  orbitals,  and  the  noded  2s  orbital  was  fitted  by  a  sum  of  three  Is  and 
three  2s  Slater  atomic  orbitals.  ^  The  fits  were  of  high  accuracy  and  were  made  for 
ease  in  the  calculation  of  the  many  integrals.  In  the  sense  that  the  energy  band  cal¬ 
culation  made  was  a  one-electron  approximation,  it  was  necessary  to  choose  an  ef¬ 
fective  potential  for  the  crystal.  The  crystal  potential  was  taken  to  be  a  superposition 
of  spherically  symmetric  atomic  potentials.  The  atomic  potential  in  turn  was  assumed 

to  be  the  Coulombic  potential  arising  from  the  effective  nuclear  charge  function,  Z  (r), 

3  p 

centered  on  the  atomic  nucleus.  The  Z„  function  for  the  P  configuration  has  been 

(3)  P 

calculated  by  Freeman,  '  and  it  was  because  of  this  convenience  that  the  particular 
configuration  of  Carbon  chosen  was  used,  since  it  was  not  felt  that  the  current  knowl¬ 
edge  of  energy  band  solutions  gave  any  other  conclusive  choice.  For  computational 

convenience,  the  Z  function  was  also  fitted,  a  linear  combination  of  four  exponentials 
P 

being  used. 

All  the  one-electron,  two-center  integrals  (overlap,  kinetic  energy  and  poten¬ 
tial)  were  done  by  the  use  of  the  usual  prolate  spheroidal  coordinate  analytic  integra- 

(4) 

tion  techniques  suitable  for  Slater  atomic  orbitals.  '  Thus,  these  integrals  required 
the  evaluation  of  auxiliary  functions^’  ^  followed  by  combinations  of  many  terms. 

All  the  one-electron,  three-center  potential  integrals  (and  also  the  two-center  poten¬ 
tial  integrals)  were  evaluated  by  the  spherical  coordinate  expansion-about-another- 

(7)  (8) 

center  technique,  similar  to  that  used  by  Lowdin,  v  '  and  Barnett  and  Coulson,  '  '  which 

(9) 

was  described  extensively  in  a  previous  Report.  ' 

It  is  perhaps  illuminating  to  give  some  of  the  computational  perspective  in¬ 
volved  in  the  overall  calculation.  Inasmuch  as  the  magnitude  of  numerical  work  in¬ 
volved  was  such  that  accuracy  would  have  been  difficult  to  maintain  in  any  hand  cal¬ 
culation,  most  of  the  computational  work  was  done  on  the  high-speed  electronic  computer 
Whirlwind  I.  The  nearly  total  mechanization  of  the  problem,  although  eliminating  al¬ 
most  entirely  any  possible  random  mistakes,  had  the  disadvantage  of  tending  to  obscure 
possible  systematic  mistakes.  The  latter  shift  of  emphasis  is  one  of  the  principle 
reasons  that  make  the  programming  of  a  computer  a  non-trivial  affair.  A  consequence 
is  that  logical  simplicity  of  the  computational  procedure  becomes  a  goal  which  is  often 
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opposition  to  computational  efficiency.  There  is  also  a  great  deal  of  difficulty  in  de¬ 
vising  adequate  test  procedures  for  computer  programs,  since  of  necessity  they  must 
be  tailored  to  the  program  itself.  The  ability  of  a  programmer  to  cope  with  these 
computer  problems  develops  mostly  with  experience.  A  large  fraction  of  the  time 
spent  on  the  present  calculation  was  thus  used  learning  how  to  obtain  the  full  poten¬ 
tiality  of  a  high-speed  computer. 

In  carrying  out  the  calculation,  the  work  fell  into  stages  for  each  of  which 
special  computer  programs  were  written.  These  were:  a  program  for  the  semi¬ 
automatic  fitting  of  the  atomic  orbitals  and  the  Z  function;  ^  a  program  for  generat¬ 
ing  the  two-center  integral  auxiliary  functions  and  then  automatically  combining  terms 
to  give  the  integrals  between  Hartree-Fock  orbitals;  ^  a  program  for  generating  the 

(9) 

atomic  orbital  expansion  functions  necessary  for  the  three -center  potential  integrals; 

a  program  for  performing  the  basic  numerical  quadratures  of  the  three-center  poten- 
(9) 

tial  integrals;  '  and  a  requantization  and  summation  program  for  forming  the  appro- 

(9) 

priate  three -center  integrals  from  the  basic  numerical  quadratures.  '  The  foregoing 
computer  programs  were  sufficient  machinery  to  prepare  the  basic  two-center  Hamil¬ 
tonian  and  overlap  integrals  which  served  as  input  for  the  final  master  program  which 
performed  the  energy  band  calculation.  Explicitly  the  basic  two-center  integrals  were 
of  the  form 


-^ij(pa)  =  /$*(?-  paL'^(r')  >pj(r )  dT 


where  the  ij>  are  atomic  orbitals,  pa  is  the  neighbor-vector,  and 


S(x) 

•M  (r)  =  i 

Jf(r) 


1  ,  (overlap) 


-jv2 


(Hamiltonian) 


(7-1) 


(7-2) 


the  sum  in  the  potential  term  being  over  all  neighbor-vectors.  In  the  present  calcula¬ 
tion,  the  three-center  potential  terms  were  neglected  when  |p^|  or  |  -  pa|  exceeded 

the  fourth  neighbor -distance  in  the  sum  over  p^. 

The  operation  of  the  master  energy  band  program  then  proceeded  as  follows. 
For  a  given  value  of  the  reduced  wave  vector,  iT,  the  program  computed  the  Hamilton¬ 
ian  and  overlap  matrix  elements  arising  from  the  Bloch  waves  constructed  from  the 
atomic  orbitals.  These  matrix  elements  were  made  real  by  taking  judicious  linear 
combinations  of  the  Bloch  waves  as  basis  states,  and  had  the  form 
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_  „  _  fcos(k-p) 

Myfk).  Z  J 

J  r-  Ism  (k  •  p  ) 

Ha 


(7-3) 


where  the  terms  in  which  p  exceeded  the  ninth  neighbor -distance  were  neglected. 

a 

The  program  next  solved  the  usual  variationally-derived  secular  equation  of  the  form 


iHj.tiT)  vjk(k)  =  E.(iT)  l  s..(iT)  vjk(k)  (7-4) 

and  stored  for  later  use  the  eigenvalues  E.(iT).  A  new  value  of  the  wave  vector  was 
then  selected  and  the  generation  and  solution  of  the  secular  equation  repeated  until 
the  pre-set  values  of  the  wave  vector  were  exhausted.  Finally  for  convenience,  the 
program  displayed  graphically  on  a  photographic  oscilloscope  cross-sectional  views 
of  the  energy  bands,  E-(iT)  vs  k,  for  values  of  the  wave  vector  along  the  edges  of  a 
basic  non-repeating  30°-60°-90°  triangle  of  the  first  Brillouin  zone. 

As  has  been  implied,  the  secular  equation  which  was  considered  had  as  basis 
states  linear  combinations  of  the  ten  one-electron  Bloch  waves  formed  from  the  Is, 

2s,  and  2p  atomic  orbitals  on  each  of  the  two  atomic  sites  in  the  spatial  unit  cell. 
Because  of  the  reflection  symmetry  in  the  Hamiltonian  operator  of  the  two-dimensional 
graphite  lattice,  the  secular  equation  one  obtains  from  these  ten  states  immediately 
factors  into  two  independent  equations,  one  of  order  eight  arising  from  reflection 
symmetric  Bloch  waves  (o-  states),  and  the  other  of  order  two  from  reflection  anti¬ 
symmetric  Bloch  waves  (ir  states).  Thus  the  master  computer  program  was  arranged 
to  independently  calculate  the  energy  band  solutions  arising  from  each  of  these  secu¬ 
lar  equations,  but  the  final  results  were  graphically  superimposed. 

The  major  physical  significance  in  a  two-dimensional  energy  band  calculation 
of  graphite  is  the  size  of  the  gap  between  the  five  lowest  (occupied)  and  the  higher 
(excited)  a  bands,  and  in  particular  whether  or  not  the  two  lower  w  bands  (valence  and 
conduction)  which  are  degenerate  at  one  value  of  the  wave  vector,  have  their  point  of 
degeneracy  in  the  <r  band  gap.  If  such  a  <r  band  gap  is  large  enough  to  include  all  of 
the  rr  bands,  then  it  follows  that  a  reasonable  approximation  for  computing  the  con¬ 
duction  properties  of  a  three-dimensional  graphite  would  be  to  ignore  the  <r  states;  in 
the  three-dimensional  crystal  the  tr  and  n  Bloch  waves  of  alternate  graphite  layers 
interact,  but  the  terminology  is  still  used.  Starting  from  the  original  work  of  Wallace, 
this  approximation  has  been  the  basis  of  all  graphite  energy  band  calculations  with  the 
exception  of  Lomer  who  has  also  considered  two-dimensional  graphite  with  <r  bands.  ^  ^ 
As  will  be  seen  later,  the  numerical  work  of  Lomer  is  believed  to  have  a  very  serious 
approximation. 
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Fig.  7-1  Fig.  7-2 

The  energy  band  solution  of  the  present  calculation  is  shown  in  Fig.  7-2 
where  the  horizontal  dimensions  are  the  edges  of  the  basic  triangle  of  the  first  Bril- 
louin  zone  as  shown  in  Fig.  7-1  .  Energy  values  of  interior  points  of  the  triangle 
were  also  found,  but  for  brevity  these  values  are  not  given  here  since  it  was  found 
that  the  interior  energy  bands  vary  smoothly  from  the  edge  band  values.  Further¬ 
more,  the  two  lowest  a  bands  arising  almost  entirely  from  the  Is  Bloch  waves, 
are  omitted  from  Fig.  7-2  since  the  bands  are  nearly  independent  of  the  wave  vector 
and  at  an  average  value  of  -15.  75  Rydbergs.  It  is  observed  that  the  highest  occupied 
point  in  the  <r  bands  occurs  at  pt.  0,  (k  =  0),  with  a  value  of  +  .  430  Rydbergs.  The 
lowest  excited  <r  band  has  roughly  a  constant  minimum  value  for  wave  vector  values 
forming  an  approximate  circle  about  pt.  0,  the  minimum  value  being  about  +  1.  158 
Rydbergs.  The  degeneracy  point  of  the  it  bands  falls  at  +  .570  Rydbergs,  a  value 
representing  the  Fermi  level  for  zero  temperature.  Thus,  it  is  clear  that  these  re¬ 
sults  support  to  some  extent  the  usual  approximation  of  neglecting  the  <r  states  in  cal¬ 
culations  of  graphite  conduction  properties.  Moreover,  this  calculation  shows  by 
virtue  of  the  overlapping  of  the  <r  and  it  bands  that  for  a  cohesive  energy  calculation, 
which  depends  on  energy  values  for  all  wave  vectors  of  the  Brillouin  zone,  the  <r  states 
must  be  included. 

Consideration  can  also  be  made  of  the  potential  function  used  in  the  present 
calculation.  The  potential,  which  was  formed  by  superposition  of  the  Zp  function, 
clearly  omits  exchange  effects.  Slater  has  given  a  procedure  for  introducing  an  ap¬ 
proximate  exchange  potential  correction  within  the  framework  of  the  one -electron 
approximation.^12^  This  exchange  potential  correction,  if  applied,  would  be  propor¬ 
tional  to  the  cube  root  of  the  charge  density  of  the  occupied  crystal  wave  functions. 
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Consequently,  one  would  expect  the  occupied  bands  to  be  lowered  more  in  energy  than 
the  unoccupied  bands,  since  the  unoccupied  wave  functions  are  orthogonal  to  the  occu¬ 
pied  wave  functions.  Thus,  it  is  plausible  that  a  more  careful  consideration  of  ex¬ 
change  effects  in  the  present  calculation  would  only  broaden  the  <r  band  energy  gap  and 
would  leave  the  present  results  qualitatively  unchanged. 

One  of  the  more  striking  features  of  the  present  results  is  the  smoothness  of 
the  energy  bands.  In  fact,  these  bands  are  similar  to  those  obtained  from  the  Slater 
and  Koster  interpolation  procedure  wherein  the  tight -binding  method  is  used  with 
neglect  of  all  but  a  few  nearest-neighbor  integrals.  It  is  therefore  of  interest  to 
examine  the  stability  of  the  present  results  in  view  of  the  several  possible  simplifying 
approximations. 

Omission  of  the  Is  Bloch  waves  was  found  to  have  an  overall  lowering  and 
warping  effect  on  the  o-  bands  in  such  a  way  that  the  o-  band  gap  was  roughly  reduced 
by  half.  Thus,  the  often-recognized  importance  of  orthogonality  is  again  emphasized. 
A  second  possible  approximation  made  was  the  omission  of  all  the  three -center  poten¬ 
tial  integrals.  This  left  the  w  bands  nearly  the  same  but  made  a  very  pronounced 
change  in  the  <r  bands,  again  closing  the  it  band  gap  down  to  about  half,  but  also  lower¬ 
ing  the  gap  so  far  that  the  tt  band  Fermi  level  no  longer  was  included. 

Finally,  the  solutions  were  examined  with  respect  to  the  approximation  of 
omitting  the  higher  neighbor  integrals.  The  stable  solution,  which  included  up  to 
ninth  neighbor  integrals,  was  found  to  be  only  slightly  warped  when  only  up  to  fifth 
neighbor  integrals  were  included  in  the  <r  bands  and  when  only  up  to  third  neighbor 
integrals  were  included  in  the  ir  bands.  However,  further  truncation  of  the  neighbor 
integrals  of  the  ir  bands  caused  violent  changes.  (It  is  for  this  reason  that  the  calcu¬ 
lation  made  by  Lomer,  which  included  only  first  neighbor  integrals,  is  not  believed 
to  be  valid.  )  In  addition,  the  effects  of  truncating  only  Hamiltonian  or  only  overlap 
neighbor  integrals  were  investigated.  It  was  found  that  the  solution  was  sensitive  to 
both  Hamiltonian  and  overlap  neighbor  truncation  to  roughly  the  same  extent  and  that 
the  two  effects  were  essentially  additive. 

Thus,  the  results  of  all  the  approximation  tests  indicate  clearly  that  the 
tight -binding  method  when  used  in  a  non-empirical  way  must  be  carried  out  with  con¬ 
siderable  mathematical  rigor  in  order  that  a  meaningful  solution  will  be  obtained. 
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8.  LATTICE  ENERGY  OF  POTASSIUM  CHLORIDE 


The  lattice  energy  of  potassium  chloride  has  been  calculated  in  conjunction 

with  a  study  of  the  electronic  structure  of  the  perfect  and  almost  perfect  crystal.  ^ 

Lattice  energy  is  the  energy  of  the  crystal  relative  to  the  energy  of  its  separated  ions; 

its  negative  is  the  cohesive  energy.  The  present  calculation  is  modeled  directly  after 

(2) 

the  comprehensive  energy  calculations  of  Dr.  P.  O.  Lowdin.  '  In  the  present  work, 
however,  interactions  between  nearest  chlorine  ions  (second  neighbors  in  the  lattice) 
have  been  included  along  with  those  between  nearest  neighbor  potassium  and  chlorine 
ions;  Lowdin  limited  himself  to  the  latter.  One  purpose  of  doing  the  new  calculation 
is  to  determine  the  contribution  to  lattice  energy  of  these  second  neighbor  interactions. 

A  further  purpose  of  the  new  calculation  is  to  provide  checks  on  quantities 
which  are  being  used  in  my  study  of  electronic  structure.  For  this  reason  my  equa¬ 
tions  for  contributions  to  lattice  energy  are  written  as  far  as  possible  in  terms  of  quan¬ 
tities  used  directly  in  the  electronic  structure  work.  Term-by-term  comparison  with 
the  calculations  of  Lowdin  and  with  calculations  based  on  the  empirical  Born-Mayer 
model  will  provide  the  check  on  these  quantities. 

Equations  for  Lattice  Energy 

The  discussion  below  follows  closely  that  given  by  Lowdin  in  his  thesis,  ^ 
and  in  general  the  notation  is  his. 

Consider  a  perfect,  infinite  crystal  of  potassium  chloride  with  its  ions  fixed 
and  in  its  ground  state.  The  crystal  is  divided  into  large  identical  regions,  a  particu¬ 
lar  region  being  symbolized  by  R.  Let  there  be  N  unit  cells  and  a  total  of  2M  electrons 
in  the  region  R.  We  take  as  a  ground-state,  many -electron  wave  function  a  single  de¬ 
terminant  $(l,  2  .  .  .  2M)  made  up  of  M  orthonormal  space  orbitals  tji.(j),  each  being 
doubly  occupied;  thus  we  take 


|(1,  2  . .  .  2M)  =  [(2M).'j  _1/2  det  |+.(j)  q(n)J  , 


(8-1) 


where  q(j)  stands  for  one  of  the  spin  functions,  a(j)  or  p(j),  and  both  occur  for  each  i. 
This  definition  leads  to  the  normalization  integral 


2M)  ${1,  l  ...  2M)  dTj  ...  dT2M  =  1  . 


(8-2) 


Use  of  the  symbol  dr  in  (8-2)  is  meant  to  imply  spin  summation  as  well  as  integration; 
use  of  dv  will  mean  integration  over  space  only. 

With  the  above  definitions  the  total  energy  of  the  crystal  becomes 


et  =  b  L  ** H  i  dTi 


c!t 


2M 


(8-3) 
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Here  H  is  the  crystal  Hamiltonian,  which  we  take  to  be  as  follows 


2M  2M 

H  =  Z  H.  +  £  G..  +  W. 

j=l  J  i.  j  =  1  3 


In  (8-4) 


V- 


’  vj 


2V  2Ze 

L  (sites  g)  — — 

1  r. 

Jg 


Gii 


M 

I* 


g.  g'  =  1 


2Z  Z  ' 
g  g 


gg' 


(8-4) 


The  letter  g  is  used  to  identify  sites,  and  is  the  nuclear  charge  on  site  g;  energies 
are  in  Rydbergs  and  distances  in  atomic  units.  Substituting  this  Hamiltonian  in  Eq. 
(8-3),  the  total  energy  becomes 


E~  = 


N 


M  M 

{ W  +  2  £  (H |  ii)  +  Z  j 2(G|  ii|  jj)  -  (G  |  ij|  ji)]}  ,  (8-5) 

i  =  1  i,  j  =  1 


where  <H  |  ii)  =  1)  Hj  vfi.(l)  dVj,  and  (G|ij|kl)  =  ^(l)  — -  4<k*{2)  i^UJdVj  dv2< 


We  now  choose  to  construct  the  orthonormal  orbitals  ijn(j)  by  making  linear 
combinations  of  the  M  normally-occupied,  free-ion,  Hartree-Fock  space  orbitals 
um(j),  centered  on  appropriate  sites  in  the  crystal.  Thus  we  write  i|r(j)  as  follows: 

M 

4*i(j)  =  Z  ami  um<3)'  3  =  1  *  • '  M-  (8-6) 

m  =  1 


The  indices  m,  n,  p  and  q  will  be  reserved  for  identifying  free-ion  orbitals.  With  this 
equation  for  t|r(j),  Eq.  (8-5)  becomes 


et  =  w 

T  N 


M  r  M  "If  M  _  M  „ 

Z  [  Z  aV  anij  (2(Himn)  +  Z  [  £  a*.a  . J[(G  |  mn|  pq)  -  (G | mq|  pn)Jj 


m,  n  =  1  i  =  1 


p,  q  =  i  j=  1 


(8-7) 


Eq.  (8-6)  relating  the  ^/'s  to  the  u's  can  be  rewritten  as  a  matrix  operation  on  vectors 
in  the  form 


¥=UA,  (8-8) 

where  ¥  and  U  are  M-dimensional  column  vectors  of  the  functions  d>.  and  u  respec- 

Ti  m 

tively,  and  A  is  the  M  x  M  matrix  of  coefficients  a  .. 

mi 

Let  us  define  an  M  x  M  overlap  matrix  A  for  the  free-ion  orbitals  by  the 
following  equation: 
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A  =  det  (Amn| ,  m  =  1  . . .  M,  n  =  1  . .  .  M,  (8-9) 


where  A 


mn 


1° 


(1)  un(l)  dVj.  Let  us  also  define 


S  =  A  -  1  ,  or  Swri  =  A  -  6 

mn  mn  mn 


(8-10) 


The  corresponding  overlap  matrix  for  the  4<'s  is  the  M  X  M  unit  matrix.  With  these 
definitions  the  orthonormality  of  the  4i's  is  expressed  by  the  matrix  equation 

l=AtAA,  (8-li) 

where  A^  is  the  adjoint  of  A.  Since  A  is  non-singular  it  has  an  inverse  A1. 
Similarly  A  has  an  inverse  A  .  Multiplying  Eq.  (8-11)  from  the  right  by  A'1  A'1 
and  from  the  left  by  A  f  (8-11)  becomes 


AA*  =  A'1  , 


or 


M 


i  =  1 


*  -  A"1 

ani  ami  “  Anm 


(8-12) 

(8-13) 


These  equations  must  hold  whether  i{>  is  a  Bloch  function  or  a  Wannier  function  or  a 
Lowdin  orthonormal  function  or  whatever.  They  follow  simply  from  the  fact  that  we 
are  constructing  M  orthonormal  functions  ij/  from  M  linearly  independent  functions  u. 

The  identity  of  Eq.  (8-13)  simplifies  the  equation  for  total  energy,  Eq.  (8-7),  and  the 
latter  now  becomes 

M  M 

Et=^+  Z  Anm  (2(H|mn)  +  Z  A‘'[2(G  |  mn|  pq)  -  (G  j  mq|  pn)]j  .  (8-14) 

m,  n  =  1  p,  q  =  1 


To  obtain  the  lattice  energy  from  E^.  we  must  subtract  the  total  energy  per 

unit  cell  of  a  crystal  with  an  infinite  interionic  distance;  the  latter  energy  we  call  E  . 

qo 

Since  the  orbitals  um(j)  are  free-ion  Hartree-Fock  orbitals  the  consistent  definition 
of  this  energy  of  the  separated  ions  is 


g 


m 


m 


Eoo=4  £  {2  ZS(Hg|mm)+  ^  ^  [2(G|mm|pp)  -  (G|mp(pm)jJ  (8-15) 

where  =  -  v2 

orbitals  u_  on  si 
m 

cept  site  g.  In  order  to  allow  explicit  subtraction  of  E  from  E^,,  we  define 

QO  1 


2Z  / r  ;  we  use  Lowdin’s  convention  that  Z  38  means  a  sum  over 
6  6  a  m 

orbitals  u  on  site  g,  while  2,  *  will  mean  a  sum  over  orbitals  on  all  sites  in  R  ex- 
m  m 


6  +  P 

nm  nm 


(8-16) 


1 
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Using  Eqs.  (8-14),  (8-15),  and  (8-16)  we  can  now  write  out  equations  for  the  crystal 
energy,  thus; 


E  -  ET  -  E«  *  Eelstat  +  Eex  +  ES’ 


(8-17) 


where 


m  p 


Eelstat  *  I  *  5  I  [>  I  '  g|g  *  I*  £*  8(G|mm| pp)j 

Eex  "sTf-  £^(G|mp|pmj),  and 

g  m  p  J 

ES  "  ^  I  I  Pnm  f2(Hl mn)  +  I  PQDf2(G|mn|pq)  -  (G|mq|pn)] 

w  g  m  n  pq  L  J 

+  2  £  [2(G|  mn  Ipp)  -  (G|mp|pn)]j. 


(8-18) 


This  separation  into  the  three  terms  of  Eqs.  (8-18)  follows  Lowdin's  procedure;  Eelstat 
is  made  up  entirely  of  electrostatic  interaction  energies,  Egx  is  a  crystal  exchange 
contribution,  and  Eg  is  a  contribution  depending  on  free -ion  overlap. 


Electrostatic  Energy 


The  electrostatic  term  in  (8-18),  Eejstat>  can  rewritten  in  a  more  obvious 
form.  Using  the  definition  of  W  from  (8-4),  we  can  write 


Jelstat 


.II  y  Z  f 

N  2  ^  fl  J 
1  g  g'fg 


2pg0)  Pg,(2) 

r  12 


dVj  dv^ 


(8-19) 


where 

p  (i)  =  Z  6(r  0)-2  ISu^(l)u  (1). 

b  &  m 


This  quantity,  Pg(I)»  is  the  charge  density  of  a  free  Hartree-Fock  ion  at  site  g,  in¬ 
cluding  its  nuclear  "charge  density".  Eeistat  *s  therefore  the  Coulomb  energy  of  in¬ 
teraction  of  Hartree-Fock  ions  superposed  to  form  a  crystal. 

For  large  interionic  distance  and  no  overlapping,  Eejstat  must  reduce  to  the 

Madelung  energy.  Let  e  be  the  valence  of  the  ion  on  site  g,  and  let  p'(l)  be  a  neutral 

g  g 

spherical  charge  distribution  of  the  ion  at  site  g,  such  that 
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PgU)  =  Pg'0)  +  *g  6(rr0)  . 


With  these  definitions, 


Eelstat  "  ^  2N^  g^gi 


Zt  t  ' 
g  g 

rgg' 


+  (j_)  i  ripg,(1)  pg,(2) 

2N  J 


g.  g' 


12 


dvldv2 


(8-20) 


The  first  term  is  directly  the  Madelung  energy,  -2a/ a;  a  is  the  Madelung  constant 
(positive  by  definition),  and  a  is  the  interionic  distance.  This  term  by  itself  accounts 
for  most  of  the  observed  cohesive  energy  and  thus  prompts  the  ionic  starting  approxi¬ 
mation.  In  the  second  term  of  Eq.  (8-20),  the  sum  over  g  can  be  reduced  to  N  times 
the  sum  over  the  two  sites  g  in  a  single  cell,  say  cell  O.  Thus  we  can  write 


E  ,  .  .  =  E  ,  +  E  , 
elstat  mad  cc 


(8-21) 


where 


and 


cc 


2 


cell  O 

Z 


g 


Z  /2pg(1)  v(2) 


dvldv2 


(8-22) 

(8-23) 


In  Eq.  (8-23)  for  E  p  (l)  and  p  ,(2)  are  neutral  spherical  distributions,  and  so  have 
c  c  g  g 

no  interaction  unless  the  electronic  shells  penetrate  one  another.  When  there  is  over¬ 
lap,  Ecc  represents  a  correction  to  Emacj  due  to  the  electronic  charge  extension;  E£c 
is  called  the  Coulomb  correction  energy  by  Lowdin.  If  the  nucleus  of  no  ion  is  appre¬ 
ciably  within  the  electronic  shells  of  any  of  its  neighbors,  then  Ecc  is  negative  and 
represents  an  attraction  of  the  ions  for  one  another. 


Exchange  Energy 

Reducing  the  sum  over  g  in  the  second  of  Eqs.  (8-18),  the  exchange  energy 


is  found  to  be 


ex 


cell  O 

I 

g 


no  g  g 

Z  Z  Z  (G|mp|pm) 


(8-24) 


m 


The  exchange  energy  is  therefore  the  sum  over  all  electrons  in  a  single  cell  of  the 
exchange  interaction  between  those  electron  and  electrons  on  other  ions.  Egx  is 
necessarily  negative  and  implies  an  attraction  of  the  ions  for  one  another. 


TheS-energy,  Ec,  depends  on  the  coefficients  P  are  derived  from  the  in- 
o  nm 
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verse  overlap  matrix,  A  * .  Lowdin  showed  that  this  matrix  A  1  can  often  be  expanded 
in  a  convergent  infinite  series  in  S (defined  by  Eq.  (8-10);  thus 

A'1  =1  -  S  +  S2  -  ...  . 

when  this  expansion  is  valid,  the  coefficients  Pnm  are  given  by 

P  =  -  S  +  V  S  S  -...  . 

nm  nm  ^  np  pm 

If  there  is  no  overlap  between  ions,  then  all  the  P  are  zero,  and  Ee  is  zero. 

nm  S 

To  give  Eg  some  meaning  let  us  look  at  the  crystal  charge  density,  p(l).  In 
our  approximation  p(l)  is  given  by 

M 

p(l)  =  2  I  +*(1)^(1)  -2  I  A^u^l)  un(l), 
i  =  1  m, n 

or  using  Eq.  (8-13)  for 

M 

p(D  =  2  I  %(l)  um(l)  +  2  £  Pnm  u^(l)  un(l)  .  (8-25) 

m  =  1  mn 

Thus  p(l)  is  made  up  of  a  superposition  of  free  ion  charge  densities,  which  of  course 
neutralizes  the  nuclear  charges,  plus  a  neutral  correction  charge  density.  The  cor¬ 
rection  term  removes  electronic  charge  from  regions  of  overlap  and  adds  it  to  the 
spherical  charge  of  the  ion.  For  KC1,  for  example,  about  0.  2  electron  charges  are 
added  to  the  original  charge  of  a  free  Cl  ion.  This  correction  charge  is  a  direct 
result  of  the  exclusion  principle  and  of  the  principle  of  the  indistinguishability  of  elec¬ 
trons;  these  principles  led  us  to  use  a  single  determinant  rather  than  a  simple  prod¬ 
uct  for  the  crystal  wave  function. 

The  S-energy  can  now  be  thought  of  as  the  self  energy  of  this  correction 
charge  density  plus  its  energy  of  interaction  with  the  charge  density  arising  from  the 
superposition  of  free  ions.  Since  Eg  includes  kinetic,  coulomb,  and  exchange  terms, 
the  contributions  to  it  are  not  all  of  one  sign.  The  overall  term  must  be  positive, 
however,  as  it  represents  a  shift  of  charge  from  the  optimum  free  ion  distribution  to 
something  more  concentrated  about  the  nucleii;  a  shift  of  this  sort  implies  that  Eg 
must  contain  terms  which  represent  an  increase  in  the  energy  of  the  "individual  ions" 
as  well  as  terms  which  reduce  the  effect  of  the  overlap-dependent,  negative  quantities, 
Ecc  and  Egx.  Both  of  these  conclusions  imply  a  positive  energy  term  Eg..  This 
S-energy  is  the  only  term  which  opposes  binding,  therefore,  and  it  is  Eg  which  pre¬ 
vents  collapse  of  the  KC1  crystal. 
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For  calculational  purposes  the  equation  for  Ec,  the  third  of  Eqs.  (8-18),  will 

S  (3) 

be  rewritten  to  involve  one -electron  matrix  elements.  From  a  previous  Report, 
these  matrix  elements  are  given  by  the  equation 

M 

Hmn  =  (Hlmn)+  £  Aqp  [^(Glmnjpq)  -  (G|mq|pn)j  •  (8-26) 

P,  q  =  1 

Using  this,  the  equation  for  Eg  (Eq.  (8-18))  can  be  written  as 

ES  =  5-  Z  I  *  I  P„m  Qmn  ■  (8-27) 

g  m  n 

where 

Qmn  =  2Hmn  '  X  Pqp  [2(GlmnIPq)  *  (G|mq|pn)]  . 
pq 

Q  is  easily  calculated  from  the  elements  H„,  and  from  terms  which  have  been 
mn  J  mn 

calculated  as  part  of  the  ingredients  of  Hmn.  In  Eq.  (8-27)  the  terms  for  which  m 
and  n  are  on  the  same  site  are  negative  while  the  others  are  positive.  An  alternative 
form  in  which  all  contributions  are  positive  is  desirable.  By  definition, 

A_1A  =!  =  (!  +  P)(l  +$)  =  1+  P+  S  +  PS,  (8-28) 


or 


P  +  S  +  PS  =  0 


The  (nm)th 


equation  of  this  form  is 


P  =  -  S  -VPS  .  (8-2*; 

nm  nm  Z.  np  pm 

P 

In  Eq.  (8-27)  for  Eg  we  can  separate  from  the  sum  those  terms  for  which  n  is  on  the 
same  site  as  m.  We  substitute  Eq.  (8-29)  for  Pnm  in  those  terms,  realizing  that  here 
Snm  is  zero;  we  interchange  orders  of  certain  summations;  we  change  the  names  of 
a  few  indices;  and  finally  we  obtain  for  the  S-energy: 


cell  O 

ES=  I 
g 


(8-30) 


In  this  form  all  the  contributions  to  Eg  are  positive  and  the  cancellation  of  terms  in 
the  square  brackets  becomes  more  exact  as  the  orbitals  involved  become  more  tightly 
bound.  Eg  will  be  calculated  from  Eq.  (8-30). 
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Calculations 

The  free-ion  orbitals  um(j)  are  chosen  real,  and  the  radial  functions  are 
linear  combinations  of  Slater  orbitals  fitted  to  the  Hartree-Hartree  numerical  tables. 
Lowdin's  fit^  was  used  for  the  important  Cl  3p  orbital,  and  new  fits  were  made  for 
the  Cl  3s  and  all  K+  orbitals.  The  orthonormalized  fitted  functions  reproduced  the 
Hartree-Fock  energy  parameters  of  Hartree  and  Hartree  to  about  the  accuracy  with 

(4) 

which  they  are  given.  '  With  these  orbitals  the  one-  and  two-electron  integrals  required 
in  Eqs.  (8-23),  (8-24),  (8-26),  and  (8-30)  were  calculated  by  the  Barnett-Coulson 
scheme  (as  adapted  for  Whirlwind  by  F.  J.  Corbato);  ^  in  this  scheme  all  the  orbitals 
in  an  integral  are  expanded  in  spherical  harmonics  about  some  one  center.  For  the 
most  part  only  integrals  involving  the  outer  shell  orbitals  (n  =  3)  were  included.  Ef¬ 
fective  nuclear  charges  for  potential  for  the  ions  were  calculated  for  use  in  Eq.  (8-23). 
Important  three-center  integrals  involving  l/r  on  one  center  were  calculated  exactly, 
while  other  three-  and  four-center  integrals  were  estimated.  The  coefficients  P 
were  found  by  using  their  symmetry  properties  and  solving  a  set  of  eighteen  simultane¬ 
ous  equations  for  eighteen  A  's  given  the  values  of  A  .  The  results  for  P  „  agreed 
with  the  first  term  of  a  Lowdin  power  series  in  S  to  one  figure;  and  they  agreed  with  a 
more  complete  expansion  to  three  or  four  figures  in  the  few  cases  tested. 

Results 

Although  E  ,  E  and  Ec  have  been  discussed  separately  above,  this  separation 

C  C  CX  d 

has  little  physical  significance.  Their  magnitudes  all  increase  with  electronic  over¬ 
lap  and  go  to  zero  if  there  is  no  overlap.  Ecc  and  Egx  are  negative,  and  the  Eg  is 
positive,  so  there  is  important  cancellation  among  them.  We  shall  therefore  lump  the 
three  terms  together  into  a  single  contribution  which  we  shall  call  the  overlap  energy, 

Eqv,  after  Born  and  Huang.  ^  To  summarize  our  equations  and  notation,  the  lattice 
energy  is  now  given  by 


E  =  E  .  +  E 

mad  ov 


where 


E  =  E  +  E  +  Ee  ; 
ov  cc  ex  S  ’ 


(8-31) 

(8-32) 


here  Ecc  is  given  by  Eq.  (8-23),  Egx  by  (8-24),  and  Eg  by  (8-18),  (8-27),  or  (8-30). 

The  results  of  the  present  calculation  are  summarized  in  row  (a)  of  the  fol¬ 
lowing  table.  The  notation  in  the  first  column  indicates  that  both  chlorine-potassium 
and  chlorine-chlorine  interactions  were  included.  The  value  of  the  interionic  distance 
a  used  was  5.  9007  a.  u.  or  3.  122  A,  appropriate  to  KC1  at  -  180  C  (where  the  V-center 
experiments  are  performed).  The  predominant  contribution  to  E  is  the  negative 
Madelung  energy;  the  overlap  energy  provides  a  ten  percent  positive  correction  to  this. 
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Table  8-1 

Crystal  energy  of  KC1:  Summary  of  results 
in  kilocalories/mole 


Source 

Emad 

E 

cc 

Eex 

ES 

Eov 

E  =  Emad 
+  Eov  U 

A)  H  (Cl,  K  +  Cl,  Cl) 

-185. 8 

-21.  0 

-60.  1 

+106. 6 

+25.  5 

-160.3  J 

B)  H  (Cl,  K) 

-185. 8 

-17.  7 

-49.  6 

(+  89.  3) 

(+22.  0) 

j 

C)  Lowdin  (Cl,  K) 

-185. 8 

- 1 6 .  5 

-48.  6 

+  84.  0 

+  18.  9 

-166.9  , 

D)  A  -  C  (Cl,  Cl) 

+  6.  5 

E)  Huggins 

-185.  8 

+  23.  2 

-162.6  -167.9  | 

F)  Exp.  (Born-Mayer)  -167.8  +  2 

It  is  simple  to  obtain  the  partial  values  of  Ecc  and  Egx,  calculated  on  Lowdin’s 
assumption  that  there  is  no  appreciable  chlorine -chlorine  overlap;  the  results  are 
shown  in  row  (B)  of  Table  8-1.  It  is  not  possible  to  find  the  value  of  Eg  appropriate  to 
Lowdin's  assumption,  however,  unless  both  Hmn  and  Pmn  are  redetermined  (Eq.  (8 -3 1)). 

It  is  possible  to  obtain  something  close  to  Lowdin's  Eg,  however,  by  extending  the  sum¬ 
mations  of  Eq.  (8-31)  only  over  pairs  (n,  m)  which  correspond  to  potassium-chlorine 
pairs.  Such  an  approximate  value  is  shown  in  parenthesis  in  row  (B)  of  Table  8-1;  a 
value  of  Eqv  calculated  from  this  approximate  Eg  is  also  shown  in  parenthesis. 

Lowdin's  results  are  shown  in  row  (C)  for  comparison  with  the  results  of  row 
(B).  The  values  of  Ecc  and  Egx  seem  to  compare  fairly  well  between  rows;  the  two 
rows  should  be  the  same,  of  course,  if  the  orbitals  used  are  identical  in  the  two  cal¬ 
culations.  Only  the  Cl  3p  orbitals  are  identical,  however,  and  the  interactions  of  my 
analytic  K+  functions  with  the  Cl  functions  might  be  expected  to  be  larger  in  magni¬ 
tude  due  to  the  longer  tails.  The  approximate  Eg  of  row  (B)  is  in  qualitative  agree¬ 
ment  with  Lowdin's  value.  The  general  agreement  between  rows  (B)  and  (C),  then, 
along  with  other  consistency  checks,  is  taken  as  providing  a  fairly  good  check  on  the 
quantities  used  here  and  to  be  used  in  the  electronic  structure  calculations. 

The  contribution  of  the  Cl,  Cl  interaction  to  lattice  energy  is  approximately 
the  difference  between  my  results  (row  A)  and  Lowdin's  (row  C);  this  would  be  exactly 
true  if  the  orbitals  used  were  the  same.  In  Lowdin's  approximation  the  overlap  energy 
is  18.  9  kcal/mole;  then  the  Cl,  Cl  interaction  contributes  another  6.  5  kcal/mole  (row 
D)  to  this  and  hence  to  the  lattice  energy.  Further  interactions  involving  inner  shell 
electrons  or  electrons  on  more  distant  neighbors  may  add  another  one  or  two  kcal/mole 
to  the  values  of  Eqv  and  E  in  row  (A),  but  probably  not  more  than  that.  Computational 
errors  probably  provide  an  additional  uncertainty  in  EQV  and  E  of  about  one  kcal/mole. 

The  value  of  E  in  row  (A)  is  therefore  probably  within  two  or  three  kcal/mole  of  the 
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correct  value  based  on  the  wave  function  assumed. 

A  comparison  of  the  foregoing  results  with  experiment  will  be  made  via  theo¬ 
retical  analysis  based  on  the  Born-Mayer  model.  The  most  up-to-date  refinements 
of  the  latter  work  were  made  by  M.  L.  Huggins.  ^  In  that  work  the  lattice  energy  of 
an  NaCl-type  lattice  was  taken  to  be  (in  C.G.S.  units) 


u  =  [-£]♦[-*  -  3]  ♦  [*KJ 


(8-33) 


x(r  +  r_  -  a) 

+  Ib6e  +  b  4.  5  e 


[» 


X(2r_  -  yia)  x(2r  -  Tia)' 

+  b  7.  5  e  r 


where  U  is  regarded  as  a  function  of  the  interionic  distance  a.  The  first  bracketed 
term  in  Eq.  (8-33)  is  the  Madelung  energy,  and  a  is  the  same  constant  which  appears 
in  (8-22).  The  second  bracketed  term  contains  the  contributions  of  van  der  Waals'  or 
dipole -dipole  interactions  between  the  ions  and  of  dipole -quadrupole  interactions;  the 
constants  C  and  D  were  determined  by  Mayer  using  crystal  spectroscopic  data.  The 
third  bracketed  term  in  Eq.  (8-33)  is  the  zero-point  energy  of  the  lattice,  v  being 
a  Debye  cutoff  frequency.  The  last  bracketed  term  contains  three  quantities,  the  first 
representing  repulsion  between  nearest  positive  and  negative  ions,  the  second  between 
nearest  negative  ions,  and  the  third  between  nearest  positive  ions.  The  constants  b, 

X,  r  ,  and  r  are  undetermined;  in  the  exponents  of  the  last  two  terms  is  the 

second  neighbor  separation;  and  the  numerical  coefficients  were  determined  by  rules 
of  Pauling. 

In  his  calculation,  Huggins  chooses  an  arbitrary  value  for  b.  He  uses  a  value 
of  X  calculated  by  Born  and  Mayer,  who  fitted  compressibility  and  equilibrium-distance 
data,  and  who  averaged  X  over  all  alkali  halides.  He  then  obtains  a  set  of  "ionic  radii" 
r+  and  r  to  predict  accurate  interionic  distances.  With  these  radii  he  recomputes  X 
and  performs  another  cycle.  The  process  leads  to  a  lattice  energy  U  which  is  the 
energy  of  the  lattice  at  0°K  relative  to  the  energy  of  the  separated  ions  at  0°K.  His 
result,  given  in  the  last  column  of  row  (E)  of  Table  8-1  agrees  with  the  experimental 
value  obtained  by  the  Born-Haber  cycle  (row  F)  to  within  0.  1  kcal/mole,  though  the 
uncertainty  in  the  latter  is  about  2  kcal/mole. 

What  has  been  called  E  in  the  present  work  is  not  the  same  thing  as  Huggins' 

U.  To  obtain  something  comparable  to  U  the  zero -point  energy  of  the  lattice  should  be 
added  to  E  and  my  interionic  distance  a  should  be  the  0°K  value.  Both  these  effects 
are  small.  The  important  difference  between  U  and  E  is  that  the  former  explicitly  in¬ 
cludes  the  dipole -dipole  and  dipole -quadrupole  interactions  between  electrons  on  differ¬ 
ent  ions  in  the  lattice.  These  effects  imply  a  correlation  of  the  motion  of  electrons  on 
different  ions.  A  single  determinant  wave  function  provides  no  possibility  for  such 
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correlation,  except  in  a  limited  way  through  exchange  correlation.  Only  by  configura¬ 
tion  interaction,  using  our  single  determinant  41  as  a  leading  term,  can  we  obtain  the 
necessary  correlation  effects. 

If  our  starting  determinant  (Eq.  (8-1))  were  the  best  possible  (determined  by 
a  Hartree-Fock  procedure),  then  separate  configuration  interactions  for  the  crystal 
and  for  the  separated  ions  should  lead  to  a  correction  to  E  which  is  roughly  given  by 
the  empirical  van  der  Waals  terms  in  Eq.  (8-33).  The  correlation  effects  within  any 
ion  should  cancel  out  in  the  process  of  subtracting  E  from  E™..  Since  our  starting 
determinant  is  not  the  best  possible,  the  configuration  interactions  must  also  correct 
for  deficiencies  in  the  one-electron  orbitals  and  yield  a  larger  correction  to  E  than  the 
empirical  van  der  Waals  terms  would  predict.  On  the  basis  of  these  remarks  it  would 
appear  that  if  our  single  determinant  is  the  best,  then  the  value  of  Eqv  obtained  from 
it  should  agree  with  the  repulsion  terms  calculated  by  Huggins  (the  last  bracketed 
terms  in  Eq.  (8-33)),  since  Huggins'  overall  result  agrees  so  well  with  experiment. 
Because  our  determinant  is  not  the  best,  Eqv  should  be  larger  than  the  Huggins  repul¬ 
sion,  leading  to  too  high  a  crystal  energy. 

The  Huggins  repulsion  terms  were  recalculated  using  his  constants,  but  my 

value  of  the  interionic  distance  a:  the  result  appears  in  row  (E)  of  Table  8-1  under  the 

heading  E  .  E  ,  in  this  row  is  of  course  the  same  as  in  all  other  rows.  The  value 
0  ov  mad 

of  Eqv  in  row  (A)  is  larger  than  the  Huggins  repulsion  by  about  ten  percent,  consistent 
with  the  discussion  of  the  preceding  paragraph.  From  our  earlier  discussion  this  dif¬ 
ference  might  increase  to  fifteen  percent  if  further  neighbors  were  included  in  the  cal¬ 
culation.  According  to  Lowdin's  results,  however,  the  calculated  value  of  Eqv  is 
substantially  less  than  the  Huggins  repulsion;  in  fact  his  value  of  E  is  almost  as  low 
as  U  itself,  and  it  leaves  no  room  for  the  correlation  effects  which  must  occur.  In¬ 
clusion  of  the  second  neighbor  Cl,  Cl  interactions  has  therefore  removed  an  apparent 
anomoly  of  the  quantum  mechanical  calculation. 

The  comparisons  above  appear  to  justify  roughly  the  starting  approximation 
of  using  a  single  determinant  of  free -ion  orbitals  for  the  crystal  wave  function.  In 
regions  midway  between  ions  this  wave  function  and  the  charge  density  it  predicts  can¬ 
not  be  too  good,  but  the  correction  charge  which  it  introduces  (Eq.  (8-25))  does  modify 
the  charge  of  superposed  free  ions  in  the  right  direction.  Small  errors  in  the  wave 
function  do  not  appreciably  affect  the  energy  value,  of  course,  and  our  value  for  the 
lattice  energy  profits  from  this  insensitivity. 

Several  interesting  facts  were  uncovered  in  analyzing  Huggins’  calculation. 

Our  equations  for  contributions  to  EQV  suggest  that  Eqv  might  vary  with  interionic 
distance  as  the  square  of  some  critical  overlap  or  overlaps.  An  overlap  integral  be¬ 
tween  orbitals  on  two  separated  ions  varies  exponentially  with  r  unless  r  is  too  small. 

A  function  A+  e  +_  r  was  therefore  fitted  to  a  curve  of  the  square  of  the  largest  Cl,  K 
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overlap  integral  versus  r,  and  a  value  for  X  of  3.  6  per  A  was  found.  Another  func- 
\  r 

tion  A__  e  "  was  fitted  for  the  largest  Cl,  Cl  overlap,  and  a  value  for  X  of  2.  8  per 

o  ~ 

A  was  found.  The  value  which  Huggins  determined  for  X  in  Eq.  (8-33),  on  the  other 

o 

hand,  was  3.  00  per  A.  Thus  Huggins'  X  is  a  mean  between  values  expected  from  vari¬ 
ation  of  these  largest  overlap  integrals. 

A  further  point  is  that  the  Huggins  value  for  Eqv>  23.  2  kcal/mole,  is  made  up 
of  three  contributions,  22.  26  for  the  +-  term,  0.  63  for  the  --  term,  and  0.  32  for  the 
++  term.  The  present  calculation  indicates  that  the  relative  contributions  of  Huggins' 

+  -  term  and  his  --  term  is  far  from  representing  the  relative  contributions  of  Cl,  K 
and  Cl,  Cl  interactions.  The  fact  that  X  is  a  mean  between  X+  and  X_  (preceding 
paragraph),  rather  than  being  much  nearer  X+  ,  again  indicates  that  the  +-  term  must 
be  including  a  large  part  of  the  Cl,  Cl  contribution.  This  fact  may  well  be  responsible 
for  the  errors  Huggins  finds  in  values  of  the  compressibility  and  of  the  maximum  op¬ 
tical  frequency  which  he  obtains  from  his  equation  for  U. 
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9.  CALCULATIONS  ON  ATOMIC  IRON 

In  writing  up  the  results  of  this  calculation,  we  have  included  the  one- 
electron  wave  functions,  one-electron  energy  parameters,  total  charge  density  for 
electrons  of  each  spin  total  potential  function,  and  the  components  of  all  the  coulomb 
and  exchange  integrals  (not,  as  yet,  kinetic  energy).  Although  it  is  expected  to  in¬ 
clude  these  quantities  in  the  published  paper,  any  and  all  are  currently  available 
from  the  writer. 

The  behavior  of  the  radii  as  a  function  of  distance  from  the  nucleus  were 

determined  from  the  charge  densities  of  each  spin,  p  and  p  .  The  size  of  these 

+  "  (l ) 

radii  are  a  measure  of  the  applicability  of  the  particular  approximation'  '  employed 
in  the  calculation. 

Let  us  consider  a  particular  orbital  u.  describing  an  electron  of  plus  spin. 
Then,  as  described  in  the  reference,  the  potential  entering  into  the  Hartree-Fock 
equation  for  this  orbital  is  constructed  from  the  total  charge  density  (p+  +  p  )  and  an 
additional  charge  distribution,  the  so-called  "Fermi”  or  "exchange  hole".  This  ex¬ 
change  hole  corrects  the  distribution  of  charge  in  the  immediate  vicinity  of  the  elec¬ 
tron  occupying  m  by  subtracting  from  p+  (in  this  case)  a  quantity  of  charge  equivalent 
to  one  electronic  charge.  The  density  of  this  quantity  of  charge  falls  off  as  we  move 
away  from  the  position  of  the  electron  and  the  calculation  of  this  density  for  the  free- 
electron  case^  shows  that  it  falls  to  small  values  a  distance  rQ  away  where 


which  would  give  the  radius  of  this  exchange  hole  if  the  exchange  charge  density  were 
uniform  and  occupied  a  sphere  just  large  enough  to  correspond  to  removal  of  one 
electronic  charge. 

One  of  the  limitations  on  the  use  of  the  free -electron  approximation  used  in 
this  calculation  for  treating  exchange  is  the  requirement  that  the  potential  in  which 
an  electron  moves  be  essentially  constant  over  the  exchange  hole.  Since  the  quanti¬ 
ties  p+  and  p  were  available,  the  radii  of  the  exchange  holes  for  plus  and  minus  spin, 
as  defined  above,  can  be  and  were  found  as  a  function  of  distance  from  the  origin. 
Considering  first  the  radius  of  the  plus  spin  exchange  hole,  we  find  that  it  has  small 
values  near  the  origin  (.  045  atomic  units’  at  .  02  a.  u.  from  the  nucleus)  and  increases 
in  a  roughly  linear  fashion  with  a  slope  of  about  1.  The  radius  for  the  minus  spin 
hole  follows  much  the  same  behavior  but  is  on  the  average  0.  3  a.  u.  greater  than  the 
corresponding  radius  of  the  plus  spin  hole. 

The  larger  radius  for  minus  spin  is  a  reflection  of  the  fact  that  we  have 
four  fewer  electrons  of  minus  spin  in  our  particular  configuration.  The  radii  are 
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somewhat  larger  than  one  might  hope  if  the  potential  were  to  be  constant  over  the 
hole. 
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10.  AN  AUGMENTED  PLANE  WAVE  METHOD  AS  APPLIED  TO 

SODIUM  AND  CHROMIUM 

Testing  of  the  matrix  element  generation  routine  is  still  going  on.  It  appears, 
however,  that  production  work  will  soon  begin.  In  the  course  of  testing  it  was  found 
that  the  normalization  constant  for  the  radial  wave  functions  for  angular  momentum 
greater  than  eight  were  negative.  The  error  was  traced  to  a  program  which,  I  believe, 
was  used  by  Howarth  in  his  calculation  on  copper.  However,  the  error  was  easily  cor¬ 
rected.  The  corrected  version  is  now  being  used  until  a  completely  new  and  recently 
written  program  for  the  calculation  of  overlap  and  ratio  of  slope  to  function  for  the 
radial  components  of  APW  is  tested. 

Work  on  chromium  has  been  halted  until  a  better  potential  than  has  been  used 
can  be  obtained.  It  is  planned  to  use  potentials  obtained  by  Nesbet's  routines  for  atoms. 


M.  M.  Saffren 


4- 


i 


11.  DESCRIPTION  OF  MANY- ELECTRON  SYSTEMS  BY  TWO-PART ICLE  ORBITALS 


I  have  described  in  previous  reports^  the  application  of  my  Coulomb  hole 
short-range  correlation  correction  to  two -electron  systems,  and  the  way  in  which  two- 
particle  orbitals  so  obtained  may  be  combined  to  form  N -electron  wave  functions.  In 
this  Report,  I  shall  mention  some  simplifications  of  the  procedure,  and  the  types  of 
integrals  which  must  be  evaluated  in  calculations  on  systems  described  by  such  wave 
functions. 

Considering  a  system  which  may  be  described  in  the  orbital  approximation  by 
a  single  determinant  $  of  orthonormal  spin-orbitals  ijj^i),  the  corrected  wave  function 
is 


$(1,  2  .  .  .  N)  =  *(1,2  .  . 


N)  (l  +  I  f’(i,j)l 
1  i<j  J 


(ll-D 


where  f'(i,  j)  are  correction  functions,  each  confined  to  the  particular  region  where 
r.j  is  small  --  the  first  restriction  being  that  regions  where  more  than  one  of  the  r^. 
are  small  are  arbitrarily  excluded.  By  the  strict  derivation  of  the  Coulomb  hole  cor¬ 
rection,  the  two-electron  function  f'(i,  j)  and  the  radius  S(i,  j)  of  the  correction  sphere 
should  depend  on  the  positions  and  spins  of  the  other  N  -  2  electrons;  the  second  re¬ 
striction  is  to  neglect  this,  and  take  f '(i,  j)  and  S(i,  j)  as  being  of  the  same  form  for  all 
pairs,  whether  singlet  or  triplet  pairs,  regardless  of  the  effect  of  the  other  N  -  2 
electrons.  It  is  by  these  two  simplifying  restrictions  that  we  can  use  the  form  ( 1 1  - 1) 
for  the  N-electron  wave  function;  it  is  also  thanks  to  these  restrictions  that  all  proper¬ 
ties  of  the  corrected  wave  function  can  be  calculated  in  terms  of  a  few  basic  integrals, 
common  to  all  systems  of  the  same  symmetry.  The  relation  of  f'(i,  j)  to  the  f(i,  j)  used 
before^  is  simply  that  f'(i,  j)  =  f(i,  j)  -  1,  so  that  f'(i,  j)  is  zero  outside  the  correction 
region.  I  further  define  G(i,  j)  =  f(i,  j)|  -  1,  which  performs  exactly  the  same  role 

in  the  expansion  of  $  or  of  as  does  ff(i,  j)  in  the  expansion  of  $  in  (ll-l),  so  that 

On  =  On{*  +  J. 

1<  j 


Using  this  expansion,  the  diagonal  density  matrices  (using  the  terminology  of 
Lowdin^)  for  the  corrected  N-electron  wave  function  can  be  expressed  in  terms  of 
the  density  matrices  for  the  uncorrected  t  and  the  following  quantities: 


^ab,  cd*1’ 
^ab,  cd^ 

^ab,  cd 


+*(1)  4*b(i)  +*<j)  G(i,  j) 

^•ij<2Sij{ab’cdll'J)  dT) 

4  .<2s.>bOi)dTi  ’  /r..<s.>, cd<i-i)dTi‘,T 
IJ  1  j  13  1J 
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the  integrals  being  taken  only  over  the  corrected  region.  The  expressions  for  the  den¬ 
sity  matrices  are  obtained  straightforwardly  by  methods  like  those  usually  used  in  ob¬ 
taining  matrix  elements  between  determinantal  wave  functions;  ^  for  instance,  the 
normalization  integral  is 

/***  drj  ...  drN  =  /*%  dTl  ...  dTN  +  2k^  ^  ^  (ll-4) 

When  the  average  value  of  a  one -electron  operator  is  to  be  computed,  two- 

electron  integrals  over  the  correction  functions  are  needed 

AS21(ab|cd)=  /  4»*(l)v|<*(2)»|id(2)  jf(l2)  Slj  *b(l)f(l2)  -0^(1)}  dr^  (11-5) 

r 12  <2S12 


For  two-electron  operators,  unfortunately  including  the  electron  repulsion  oper¬ 

ator  in  the  Hamiltonian,  three-electron  correction  integrals  are  needed: 


An 


12 


(ab  |  cd  |  ef)  =  J  +*(l)  **(2)  +*(3)  +f(3)  {f(l3)012  +b(l)  ^(2)  f(l3) 

r13<2S13 

-  «124>b(l)  +d(2)}  dr1dT2dT3 


(11-6) 


I  am  now  writing  programs  to  compute  the  one-  and  two-electron  integrals 
like  ( 1 1  - 3 )  and  ( 1 1  -5)  for  the  case  of  spherical  symmetry  --  that  is,  for  the  He  atom, 
the  same  integrals  with  changed  parameters  being  applicable  to  any  atom.  The  case 
of  atoms  not  in  S  states  can  be  handled  by  the  same  programs  by  the  device  of  intro¬ 
ducing  an  angular  part  of  the  correction  function  f(i,  j)  so  that  the  total  corrected  wave 
function  in  the  corrected  region  is  made  independent  of  the  direction  of  the  interelec¬ 
tron  vector  --  the  essential  dependence  on  the  magnitude  r.j  thus  not  being  affected. 

The  three -electron  integrals  like  (11-6)  are  not  yet  being  attempted.'  Com¬ 
putational  difficulty  for  a  one -center,  three -electron  integral  is  comparable  with  that 
for  a  two-center,  two-electron  integral. 
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12.  ON  THE  POLARON  ENERGY  AND  MASS 

R.  P.  Feynman^  ^  in  his  approach  to  the  polaron  problem  has  compared  the 
problem  of  an  electron  interacting  with  the  phonon  field  with  that  of  an  electron  inter¬ 
acting  harmonically  with  a  second  particle  of  arbitrary  mass.  He  has  achieved  a 
variational  expression  for  the  ground  state  energy  of  the  first  problem  in  terms  of  the 
two  variable  parameters  of  the  second,  the  mass  M  of  the  second  particle  and  the 
spring  constant  k  of  the  harmonic  coupling.  Rather  than  these  two  parameters,  Feyn¬ 
man  introduces  the  parameters  v  and  w  related  to  M  and  k  by 


=  Vk/M  and  v  =  ->/k(l  +  l/M) 


(12-1) 


In  terms  of  these  parameters,  an  upper  bound  to  the  energy  is  given  by 


EF<Vw>  "  *F<V  -  w)2  -  / 


e  *  dt 


^  Jo  rj~ 

V  w  t  + 


(units  of  hu)  (12-2) 


2  2 

v  -  W  /,  -vt, 
- (1  -  e  ) 


which  is  to  be  minimized  with  respect  to  v  and  w.  For  the  best  choice  of  v  and  w  an 
estimate  of  the  effective  polaron  mass  for  very  slow  polarons  is  given  by 


mF  =  1  + 


av 


r 


-t  ,2  .. 
e  t  dt 


3y^  Jo  /  2  v2  _ 

v  w  t  +  - 


(12-3) 


w 


■0  -  e-vl) 


Feynman  obtains  the  latter  expression  by  making  what  amounts  to  a  kind  of  first  order 

correction  to  the  total  mass  of  the  approximating  two  particle  system,  m~  =  M  +  1  = 

2  2  1 
v  /w  .  Feynman  also  has  given  approximate  expansions  for  v,  w,  Ep,  and  mF  in  the 

weak  and  strong  coupling  limits.  Because  of  the  conflicting  estimates^  of  the  polaron 
mass  in  the  coupling  range  a  =  5  -  7,  it  was  felt  advisable  to  carry  out  the  minimiza¬ 
tion  of  (l2-l)  and  the  determination  of  the  mass  to  zero  order,  m^,  and  to  first  order, 
m,,,  comparing  these  with  each  other  and  with  masses  calculated  by  other  methods. 

The  repeated  evaluation  of  ( 12- 1)  for  many  values  of  the  pair  (v,  w)  was  per¬ 
formed  by  the  Whirlwind  computer.  Any  program  for  the  evaluation  of  these  integrals 
had  to  contend  with  the  singular  behavior  of  the  denominator  for  t  -*  0,  the  rapidly 
changing  factor  e  vt  in  the  denominator,  and  the  infinite  range  of  integration.  The 
program  finally  developed  treated  the  range  in  two  parts,  from  zero  to  t  and  from 
t  to  co,  where  t  *  7.  5/v.  In  the  first  region  the  variable  transformation  y2  =  t/tQ 
removes  the  singularity  at  the  origin,  includes  all  the  range  of  troublesome  behavior 
of  e~vt  in  a  way  independent  of  v,  and  allows  accurate  integration  with  16-point  Gauss- 

Legendre  quadratures.  ^  In  the  second  region  the  variable  transformation  y  =  t  -  t 

(4)  ° 

allows  accurate  integration  with  15-point  Gauss -Leguer re  quadratures. '  '  Since  the 
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energy  is  very  insensitive  to  the  choice  of  v  and  w,  six-place  accuracy  in  computing 
( 1 2  - 1)  was  required.  This  was  achieved  in  a  computer  code  written  by  Miss  Hannah 
Paul  requiring  between  one  and  two  seconds  per  integral.  Subsequent  minimization  of 
E(v,  w)  was  carried  out  graphically  by  Miss  Paul. 

In  Table  12-1  we  have  listed  the  values  of  v,  w,  Ep,  m^,,  and  mp  found  for 
various  values  of  a.  For  comparison  we  have  also  listed  the  energy  (E.  T  p  =  -  o)  and 


mass  (m 


LLP 


(51  LLP  j  + 

=  1  +  a/6)  found  by  Lee,  Low  and  Pines  and  Gurari,  '  '  the  mass(mUt>-  ; 


(6) 


Ha  i 


found  by  Haga'  '  based  on  the  LLP  procedure,  the  energy  we  have  computed  using 
Gross's  transformation  superposed  on  the  LLP  shift  function,  ^  and  the  energy  (Epp>T 
-  .  1088a2  -  3/2)  and  mass  (m„p„  =  232(a/l0)4)  obtained  by  Pekar,  ^  Bogolyubov, 

(  1  Q\ 

and  Tyablikov.  '  It  should  be  noted  that  the  Pekar  energy  expression  is  not  derived 
from  a  variational  calculation  so  that  it  may  be  lower  than  the  true  energy.  It  is  true 
that  Tyablikov^1  ^  and  Hohler  have  derived  this  expression  variationally  but  only  in  the 
limit  of  very  large  coupling. 

What  is  immediately  noticed  is  that  the  Feynman  energy  is  significantly  lower 
(hence  better)  in  the  coupling  region  of  interest  than  the  energy  from  any  other  theory 
and  that  the  Feynman  masses  MT  and  Mp  differ  considerable  from  both  the  LLP  and 
PTB  values  in  this  range  of  a.  Furthermore,  the  behavior  of  Haga's  expression  for 
a  *♦  12  is  again  seen  to  be  objectionable. 

In  addition  to  the  fact,  shown  by  Feynman,  that  m„  has  essentially  the  right 

*  * 

behavior  in  both  the  weak  and  strong  coupling  limits,  we  now  have  the  encouraging 
information  that  mp  makes  only  a  small  correction  to  the  simpler  m,p  at  all  coupling 
strengths.  Thus,  a  theory  of  mobility,  for  example,  which  contains  implicitly  the 
approximation  of  m^,  for  the  polaron  mass  should  not  be  so  much  in  error.  Theories 
using  mLLp  or  nipgrp  may  be  considerably  in  error  if  only  because  of  an  incorrect 
description  of  the  slow  polaron. 


References 

1.  R.  P.  Feynman,  Phys.  Rev.  97,  660  (1955). 

2.  T.  D.  Schultz,  Quarterly  Progress  Report,  Solid-State  and  Molecular  Theory 

Group,  M.  I.  T. ,  October  15,  1955,  p.  43. 

3.  A.  N.  Lowan,  N.  Davids,  A.  Levinson,  Bull.  Am.  Math.  Soc.  48,  739  ( 1 942). 

4.  H.  E.  Salzer  and  R.  Zucker,  Bull.  Am.  Math.  Soc.  55,  1004  (1949). 

5.  T-D.  Lee,  F.  E.  Low,  and  D.  Pines,  Phys.  Rev.  90,  297  (1953);  M.  Gurari, 

Phil.  Mag.  44,  329  (1953). 

6.  E.  Haga,  Prog.  Theor.  Phys.  (Japan)  13,  555  (1955). 


$ 

The  small  constant  percentage  error  occurring  as  a  -*■  ao  is  due  essentially  to.the 
restricted  shape  of  the  effective  potential  well  required  by  the  harmonic 
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13.  ON  THE  MOBILITY  OF  THE  POLARON 

The  simplest  property  of  the  polaron  which  one  has  attempted  to  observe  ex¬ 
perimentally  is  not  the  self  energy  or  effective  mass  of  a  slow  polaron  but  the  polaron 
mobility.  Attempts  have  been  made  to  calculate  the  electron  mobility  in  polar  crystals 
including  some  treatment  of  self  energy  and  mass  corrections  by  Frohlich,  Pelzer  and 
Zienau,  Low  and  Pines,  ^  Morita  et  al.  ,  ^  Pekar^  and  the  author.  ^  Character¬ 
istic  of  these  approaches  has  been  the  attempt  to  base  the  mobility  calculation  on 
standard  transport  theory  applied  to  a  more  subtle  description  of  the  polaron  than 
calling  it  a  bare  electron.  More  direct  attempts  to  calculate  the  mobility  without  a 
detailed  calculation  of  free  polaron  properties  or  a  specific  reliance  on  the  Boltzmann 
equation  have  been  made  by  Lax^  and  Feynman,  ^  although  these  have  not  yet  pro¬ 
vided  results  comparable  to  the  earlier  works.  Of  these  more  complete  analyses, 
only  that  of  Morita  et  al.  does  not  rest  on  one  particular  description  of  the  polaron, 
and  those  of  FPZ,  LP,  and  Pekar  depend  on  descriptions  which  seem  inferior  to  the  Feyn¬ 
man  approach  (see  preceding  Report).  It  has  therefore  been  our  aim  to  apply  the  Feynman 
description  of  the  polaron  to  the  mobility  problem.  An  approximate  application  of 
this  technique  has  been  achieved  and  will  be  described  in  detail  in  the  author's  doc¬ 
toral  dissertation  and  eventually  in  the  literature.  We  shall  restrict  ourselves  here 
to  summarizing  a  few  salient  features. 

(1)  The  phonon-polaron  scattering  problem  is  best  considered  as  a  resonance 
scattering  problem  for  two  reasons.  First,  it  avoids  the  necessity  for  a  description  of 
the  motion  of  a  "fast"  polaron,  i.  e.  one  which  has  absorbed  a  phonon  and  so  can  emit 
one  spontaneously.  Such  a  polaron  is  not  even  approximately  in  a  stationary  state  in 
most  cases,  and  cannot  be  described  in  the  same  way  as  a  polaron  with  small  but  non¬ 
zero  momentum.  Second,  it  avoids  the  difficulties  discussed  by  Howarth  and  Sond- 

(8) 

heimer'  '  in  solving  the  Boltzmann  equation  when  the  simple  scattering  process  pro¬ 
duces  a  significant  change  in  the  particle  energy. 

(2)  It  may  well  be  that  the  resonance  scattering  can  be  considered  as  that  of 
a  field  quantum  by  a  quasi-particle,  but  even  if  this  is  so  one  must  know  not  only  the 
energy-momentum  relation  for  the  particle  near  zero  momentum,  but  also  the  reson¬ 
ance  momentum  and  velocity  (calculated  from  knowing  the  energy -momentum  relation 
well  away  from  zero  momentum)  and  the  width  of  the  intermediate  (resonance)  state. 

In  this  respect  the  assumptions  of  FPZ  and  of  Pekar  prove  unacceptable. 

(3)  There  is  evidence  from  LP's  work  that  the  resonance  scattering  cannot 
be  considered  simply  as  that  of  a  field  quantum  by  a  quasi-particle  as  manifested  in 
the  (m/m  )  rather  than  (m/m  )  '  dependence  of  mobility  on  polaron  mass  found  in 

their  paper.  The  analysis  of  Morita  et  al.  however  suggests  to  the  contrary. 

(5) 

(4)  The  simple  replacement,  described  in  a  previous  Report,  '  '  of  the  exact 
effective  action  for  the  electron  propagation  by  the  approximate  action  introduced  by 
Feynman  is  unsatisfactory  in  one  important  respect.  Since  the  effective  action  describes 


St  . 
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the  motion  of  one  particle  in  a  two-particle  system  having  only  bound  states,  it  is  in¬ 
capable  of  giving  a  finite  width  to  the  polaron  state  after  absorption  of  a  phonon.  Since 
this  width  is  crucial  in  obtaining  a  finite  resonance  scattering  probability,  the  results 
previously  reported  are,  as  they  stand,  meaningless.  Furthermore,  any  finite  order 
correction  to  the  use  of  the  simple  approximate  action  obtained  by  summing  over  a 
finite  number  of  Feynman  diagrams  with  virtual  phonons  representing  the  difference 
between  the  approximate  and  exact  effective  actions  cannot  give  a  finite  width  to  this 
intermediate  state  just  as  no  finite  order  of  perturbation  theory  could  give  the  Weisskopf 
resonance  fluorescence  formula  in  light  scattering  by  atoms. 

(5)  A  special  but  infinite  set  of  Feynman  diagrams  can  be  summed  over  in  this 
case  with  good  approximation  and  a  width  is  obtained.  These  calculations  are  best  per¬ 
formed  by  proving  the  equivalence  of  the  Feynman  formulation  of  the  polaron  problem 
with  an  extended  Hamiltonian  formulation  and  then  doing  the  sums  in  the  Hamiltonian 
formulation.  One  finds  that  the  picture  of  resonance  scattering  as  that  of  a  field  quan¬ 
tum  by  a  quasi -particle  is  valid,  at  least  to  this  approximation,  contrary  to  the  results 
of  Low  and  Pines,  and  in  agreement  with  Morita  et  al.  That  is,  the  mobility  depends 
on  the  polaron  mass  as  (m/m  )  '  .  One  also  finds  that  if  similar  corrections  are 
calculated  for  the  slow  motion  of  a  free  polaron,  the  effective  mass  is  changed  by  only 
a  few  percent  from  m^,  or  mf  discussed  in  the  previous  Report. 
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